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Summary 

The  project  addressed  a  number  of  fundamental  aspects  electromagnetic  fields  and  their 
applications  in  inverse  scattering  problems.  The  program  advanced  the  understanding  of 
generation,  propagation,  manipulation  and  detection  of  random  electromagnetic  beams. 
The  research  effort  comprised  theoretical,  numerical  modeling,  and  experimental 
developments. 

In  order  to  simulate  the  results  of  scattering  from  a  random  medium  to  a  degree  of 
accuracy  that  will  reflect  the  statistical  nature  of  the  studied  phenomenon,  the  coupled 
dipoles  approximation  (CDA)  is  a  method  of  choice.  CDA  does  not  suffer  from  the 
important  disadvantage  of  FDTD  and  FEM,  namely  the  need  to  discretise  space  outside 
the  medium  and  the  need  to  implement  suitable  boundary  conditions  to  prevent  non¬ 
physical  reflections  from  the  boundaries  of  the  computational  domain.  We  have 
developed  a  CDA  algorithm  capable  of  dealing  with  different  realizations  of  a  random 
medium  and  statistical  analysis  of  the  light  scattered  as  a  result  of  light  matter  interaction. 

In  most  cases  of  interest  for  remote  sensing,  the  result  of  interaction  between  coherent 
waves  and  random  media  is  scale  dependent.  We  investigated  the  coherent  scattering 
from  random  media  both  numerically  and  experimentally.  We  have  developed  a 
formalism  to  describe  the  vectorial  properties  of  random  media  and,  based  on  the  scale 
dependent  responses,  we  have  suggested  unique  sensing  approaches.  Our  findings  are 
relevant  for  both  material  identification  and  the  design  of  novel  materials  because  we 
have  identified  a  new  electromagnetic  interaction,  which  represents  the  scale  at  which  the 
polarimetric  response  of  a  medium  is  most  sensitive  to  the  excitation  field 

We  have  demonstrated  a  new  way  of  analyzing  the  fluctuations  of  scattered  waves  and 
have  shown  that  individual  members  of  an  ensemble  of  interactions  provides  means  to 
extract  information  beyond  that  available  in  the  ensemble  average.  This  should  be  of 
particular  interest  for  remote  sensing  procedures  that  rely  on  multiple  exposures  or 
extended  time-scales  in  order  to  establish  robust  averages. 


We  have  also  demonstrated  that,  in  spite  of  usually  being  considered  a  nuisance,  media 
causing  intricate  scattering  can  act  as  efficient  linear  optical  devices.  We  have  shown 
that,  upon  calibration,  random  media  can  be  used  as  polarimeters  with  remarkable 
properties. 

When  coherent  beams  interact  with  random  media,  a  number  of  polarization  dependent 
phenomena  occur.  For  instance,  due  to  spin-orbit  interaction,  the  optical  wave  encounter 
with  a  gradient  of  refractive  index  leads  to  transport  of  spin  similar  to  the  electronic  spin 
Hall  effect.  We  have  demonstrated  that  for  a  wave  in  a  pure  state  of  polarization,  the 
spin-orbit  interaction  results  in  a  spiraling  power  flow  that  is  determined  by  the  extent  of 
interaction.  This  constitutes  the  first  demonstration  that  spin  transport  can  be  manipulated 
and  enhanced  in  confined  geometries.  The  results  open  the  possibility  for  new 
functionalities  for  sensing  devices  and  unique  communication  channels. 

It  has  been  known  for  quite  some  time  that  electromagnetic  fields  can  induce 
conservative  forces  resulting  from  field  gradients  as  well  as  nonconservative  forces 
appearing  due  to  radiation  pressure  and  gradients  of  phase.  We  have  uncovered  a  novel 
interplay  between  conservative  and  nonconservative  forces,  which  constitutes  a  new 
mechanism  to  induce  torques  on  spherically  symmetric,  optically  isotropic,  and  lossless 
objects.  The  phenomenon  can  be  exploited  to  generate  torques  that  are  fully  controlled  by 
the  polarization  of  the  excitation  field.  Because  the  torques  acting  on  components  of  a 
material  system  are  entirely  nonconservative,  the  phenomenon  represents  a  new  means  to 
extract  energy  form  a  radiation  field  besides  the  traditional  dissipation  mechanisms. 

Certain  results  of  this  program  have  been  included  in  a  number  of  publications:  two  Phys 
Rev  Lett,  two  Phys  Rev  E,  and  three  Opt.  Exp.  The  program  has  also  supported  a  number 
of  twenty  three  presentations  and  invited  presentations  at  national  and  international 
scientific  meetings. 


l.  Coupled  dipole  approximation  for  random  media 

The  coupled  dipoles  approximation  (CDA)  is  a  numerical  technique  that  discretizes  a 
continuum  volume  in  a  finite  array  of  polarizable  point  dipoles  [1].  These  dipoles  react 
to  the  local  field  and,  by  the  use  of  dyadic  Green's  function,  the  complete  interactions 
between  dipoles,  including  both  the  near-  and  far-field  components,  are  accounted  for. 
The  polarizability  of  each  individual  dipole  relates  to  the  local  dielectric  properties  of  the 
medium.  In  modeling  a  near  field  experiment,  the  initial  excitation  is  produced  by  an 
emitting  dipole  placed  in  the  near  field  of  the  medium  which  simulates  the  tip  of  the 
fiber. 

According  to  coupled  dipole  approximation  the  field  at  each  of  the  dipoles  can  be 
written  as  the  summation  of  both  the  incident  field  and  the  contributions  occurring  from 
the  interaction  with  all  the  other  dipoles  on  the  lattice 

E(rj)  =  Einc(rj)  +  fjEkdip(rj),  (1) 

k=  1 

where  Ekd,p(rj )  is  the  electric  field  radiated  to  the  point  r  by  the  kth  point  dipole  Pk . 
Considering  a  point  dipole  at  an  arbitrary  orientation  and  at  an  arbitrary  location  rk ,  the 
radiated  field  can  be  written  as: 

rt)\-Pt,  (2) 

where  A(r.,^)is  the  interaction  matrix  accounting  for  both  near-  and  far-field  coupling 

components  as  a  result  of  the  dyadic  Green's  function.  This  matrix  is  a  dense,  symmetric 
positive  definite  matrix  of  size  3Nx3N  where  N  is  the  number  of  dipoles  used.  The 
elements  of  A{rj,rk)  can  be  found  by: 
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where  rjk  =  r-  rk 


and  rjk  = 


To  model  inhomogeneous  media  using  CDA,  dipoles  with  polarizabilities  given 
by  the  Lorenz-Lorenz  formula  corresponding  to  the  specific  components  were  randomly 
distributed  in  a  proportion  according  to  the  sample's  composition.  The  far  field  coherent 
superposition  of  the  dipolar  contributions  is  recorded  for  specific  configurations  of  the 
scattering  medium  and  an  ensemble  of  realizations  is  produced.  A  large  number  of 
realizations  of  a  sample  are  generated  and  the  relevant  statistical  information  of  phase, 
amplitude,  and  intensity  distributions  are  determined. 


1.2  Random  Materials  and  CDA 

In  the  specific  case  of  modeling  a  random  composite  material,  dipoles  with  the 
different  polarizabilities  are  randomly  distributed  across  the  lattice  in  a  proportion 
determined  by  the  medium's  composition.  In  our  procedure,  this  random  distribution  of 
polarizabilites  constitutes  one  realization  of  the  randomly  inhomogeneous  medium.  An 
ensemble  of  such  realizations  is  obtained  by  creating  a  large  number  of  independent 
realizations.  Figure  shows  an  example  of  spheres  randomly  distributed  throughout  the 
cubic  lattice. 
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Figure  1  One  realization  of  spheres  imbedded  in  a  host  material.  A  sphere  on  a  cubic 
lattice  in  this  situation  is  modeled  by  7  dipoles  and  is  randomly  distributed  throughout  the 
lattice  according  to  the  volume  fraction. 

The  number  of  iterations  needed  to  approach  the  solution  depends  on  both  the 
refractive  index  contrast  between  the  different  constituents,  and  on  the  initial  guess  at  the 
solution.  For  each  realization  of  the  sample,  the  far-field  coherent  superposition  of  the 
dipolar  contributions  is  recorded.  After  a  large  number  of  these  realizations,  an  intensity 
distribution  is  created  and  statisticaly  relevant  information  such  as  the  optical  contrast  can 
be  analyzed. 


Figure  2.  Distribution  of  states  of  polarization  for  field  scattered  from 
simulated  cube  of  random  medium.  Blue  points  show  experimental 
constraints,  red  dots  show  real  Stokes. 


1.3  Extension  of  CD  A  to  infinite  media 


One  of  disadvantages  of  coupled  dipole  approximation  is  the  rapid  increase  of 
memory  requirements  with  increase  of  the  number  of  dipoles.  One  of  the  possible  ways 
to  avoid  this  problem  is  to  introduce  periodic  conditions  for  modeling  geometry.  For 
example,  to  study  the  scattering  properties  of  infinite  slab  of  inhomogeneous  material 
using  CDA,  the  modeling  cube  can  be  replicated  in  2  dimensions,  (x-y  plane)  with  some 
periodicity  specified  by  the  initial  cubes  with  a  side  dimensions  d.  The  assumption  of  this 
extension  is  that  this  periodicity  will  not  modify  the  essential  properties  of  the  scattered 
light. 

We  would  like  to  consider  the  situation  where  we  have  a  plane  wave  with  some 
arbitrary  angle  of  incidence  on  the  slab  of  inhomogeneous  material: 
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Emc  =  Eo(r.)e 


(4) 


According  to  the  periodicity  condition  for  such  an  excitation,  we  following  relation 
should  be  satisfied  for  the  field  inside  the  slab: 


E(rj  +  ndx  +  mdy)- 


(5) 


Here  x  and  y  are  unit  vectors  in  the  x-  and  y-directions,  kt  is  the  component  of 


wavevector  parallel  to  the  surface  of  the  slab.  The  initial  equation  for  coupled  dipole 
approximation  (1)  now  can  be  rewritten  as  follows: 

E{Vj)  =  E0(fj)+  Yj  ajEirj)  YA(rj,rk  +  ndx  +  mdy)e,k"{ndx+mdy) 

rteVd  m,n=-cc  ^ 

The  first  summation  in  the  right  hand  side  of  this  expression  is  performed  over  all  dipoles 
in  the  modeling  cube.  The  second  summation  is  performed  over  each  monolayer 
corresponding  to  the  repetition  of  each  plane.  Unfortunately,  the  lattice  sums  over  n  and 
m  in  real  space  does  not  converge.  However,  these  lattice  sums  may  be  calculated  after 
applying  two-dimensional  Fourier  transformation.  After  proper  transformations  we  have 
the  following  expressions  for  the  lattice  sums: 
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In  this  case  the  sum  in  real  space  does  not  converge,  and  the  sum  in  Fourier  space 
diverges,  as  such,  it  is  necessary  to  use  another  strategy.  In  this  case,  the  sum  can  be 
written  as  a  combination  of  the  sum  in  real  space  and  reciprocal  space.  The  convergence 
is  made  faster  by  adding  and  subtracting  and  offset  sum  as  a  function  of  some  arbitrary  z 
value  of  h  from  the  real  space  sum 


y  A  (r. ,  rk  +  ndx + mdy)e'k^ri+ndi+mdy^  =  a(kt, hj+  {a  (ijj ,o)—  A (ij , hjj 
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Sum  in  inverse  space  is  the  same  as  the  first  case  however  now  the  z  separation  is  offset 
from  zero  to  be  some  value  h: 
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And  the  sum  in  real  space  is  the  usual  dipole  interaction  Green  function,  evaluated  at  z=h 
and  z=0 
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0  in  formula  (9)  denotes  the  outer  product, 

Kq  =  +  8Pm  ’  si8n{zj  ~  zk  \Pq\  gm  (11) 

*pq  =  Jk*-K+lJ  ’  Pu=fc-?j\ 

where,  p  and  q  are  integers. 

These  new  expressions  should  allow  the  possibility  to  model  a  semi-infinite  slab 
without  the  influence  on  the  boundary  effects.  However,  the  new  expressions  are  no 
longer  symmetric  with  respect  to  the  separation  vector.  As  CDA  involves  solving  a  large 
dense  system  of  equations,  it  often  necessary  to  use  numerical  techniques  that  may 
drastically  improve  the  efficiency  of  obtaining  this  solution.  The  most  common 
technique  is  the  use  of  a  Fourier  Transform  technique  and  an  iterative  approach  to  the 
solution  known  as  the  conjugate  gradient  algorithm.  This  technique  draws  on  the 
symmetries  available,  and  as  such,  our  new  Greens  function  should  be  symmetric. 

The  system  can  become  symmetric  by  multiplying  both  sides  by  e~,k  p  and 
solving  a  new  system  of  equations 

=  E0{rj)e~lk'p  +  ^  ^  A  (r ,  rk  +  ndx  +  mdy]e'k^"dx+mdi'ajE(rj )e~lk'p 

rj<EVd  m,n=-oo  (I'l) 

One  feature  that  has  not  been  discussed  thus  far  is  the  influence  of  the  periodicity  on  the 
scattered  fields.  In  general,  the  scattering  of  a  plane  wave  from  a  periodic  structure 
results  in  diffracted  orders  that  depend  on  the  wavelength  and  the  periodicity  length. 
Figure  3  shows  directions  of  diffracted  orders  for  different  periodicity  parameters  for 
incident  field  at  45°. 


Figure  3.  Directions  of  diffractive  orders  for  periodicity  parameter 
equal  to  d  =1.6 A,  (left),  d  =3.2/1  (right). 

It  is  possible  to  show  that  in  the  case  of  homogenous  material,  only  the  0th  order 
diffraction  occurs,  corresponding  to  the  reflected  wave.  In  the  case  of  inhomogeneous 
material,  these  diffracted  orders  do  not  disappear.  But  the  presence  of  additional 
diffractive  orders  is  not  the  artifact  of  method.  Let  us  suppose  that  the  periodicity 
approaches  infinity.  Than  direction  of  diffractive  orders  occupy  the  whole  space.  Thus, 
this  additional  scattered  light  that  appears  in  simulations  could  serve  as  some  sort  of 
measure  of  diffuse  light. 

Initial  testing  of  the  method  was  done  in  situation  where  the  material  is 
homogenous.  For  this  situation,  there  exist  analytical  solutions  for  the  field  inside  the 
slab  as  well  as  the  corresponding  reflection  coefficients.  In  Figure  4  we  present  the 
comparison  between  the  field  distribution  inside  the  slab  calculated  analytically 
according  to  Fresnel  formulas  and  the  results  obtained  using  our  CDA  extension.  The 
calculations  demonstrate  an  excellent  correspondence  between  the  two  approaches. 
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Figure  4.  Electric  field  magnitude  inside  plane-parallel  slab  (64  dipoles 
in  depth,  total  0.8X).  Blue  curve  shows  analytical  results,  red  dots  show 
the  result  obtained  with  CDA  extension  to  infinite  media. 


The  new  method  was  further  developed  to  deal  with  modeling  large  scale  random  media. 
We  have  demonstrated  that  this  method  can  adequately  describe  optical  properties  such 
as  specular  and  diffuse  reflection  amplitudes.  The  accuracy  of  our  approach  was  tested  on 
slabs  with  thicknesses  less  than  the  wavelength  of  light.  We  found  an  excellent  agreement 
with  analytical  results  for  slabs  of  homogeneous  media  and  for  media  containing  small 
size  inclusions.  When  the  size  of  inhomogeneities  increases,  our  numerical  approach  can 
describe  phenomena  for  which  no  analytical  models  are  available. 


Z 


Fig.  5  Infinite  slab  of  an  inhomogeneous  material  modeled  using  CDA.  The 
boundary  effects  are  eliminated  by  replicating  a  random  volume  in  2D. 

Further  expansion  could  also  address  the  problem  of  semi-infinite  inhomogeneous 
media.  However,  in  this  case  a  direct  extension  of  our  approach  is  not  possible,  because 
the  electromagnetic  field  in  z-direction  is  not  exactly  periodic.  The  presence  of  the 
surface  breaks  the  symmetry  and  makes  the  field  in  the  surface  region  different  from  the 
field  in  the  bulk.  The  generalization  procedure  should  take  this  into  account. 

Another  particular  direction  of  interest  is  modeling  roughness  effects  by  introducing 
irregularities  into  the  upper  surface  of  the  modeling  cube.  In  the  case  of  coupled-dipole 
method,  media  with  irregular  surface  structure  represent  a  particular  case  of 
inhomogeneities  where  some  surface  dipoles  have  zero  polarizability.  Thus,  the 
formalism  described  above  for  inhomogeneous  media  is  valid  also  for  media  with  rough 
surfaces. 


2.  Scale-dependent  description  of  inhomogeneous  media 

A  second  effort  was  directed  at  finding  means  to  describe  the  effective  polarizability  of 
random  media.  When  describing  material  properties,  the  probed  scale  determines  the 
outcome.  The  observed  material  response  at  certain  scales  can  depend  on  the  local 
structure.  In  general,  one  can  define  microscopic,  mesoscopic,  and  macroscopic  regimes, 
but  their  absolute  length  scales  will  depend  on  the  specific  material  characteristics. 

The  microscopic  scale  refers  to  the  smallest  volume  over  which  the  material  is 
homogeneous;  probing  this  scale  reveals  the  “intrinsic”  properties  of  the  medium.  The 
optical  characteristics  at  the  largest,  macroscopic  scale  involve  significant  averaging  over 
volumes  usually  much  larger  than  characteristic  scales  of  inhomogeneity  in  medium. 
Mesoscopic  scales  on  the  other  hand,  are  small  enough  that  fluctuations  around  the 
average  become  important  and  may  contain  additional  information  about  structural 
characteristics  of  inhomogeneity.  In  this  regime,  the  classical  laws  of  interaction  do  not 
apply  anymore. 

A  material  is  optically  homogeneous  when  its  optical  responses  are  equivalent 
from  micro-  up  to  macroscopic  scales.  In  the  case  of  so-called  inhomogeneous  media,  a 
microscopic  description  provides  the  intrinsic  dielectric  properties  of  constituents;  while 
a  description  at  macroscopic  scales  results  in  an  effective  averaging  of  the  dielectric 
properties.  All  known  effective  medium  theories  are  based  on  this  averaging  principle1. 


Figure  6.  (a)  A  material  system  with  varying  optical  properties  over  some 

characteristic  mesoscopic  length  scale.  (b)  Conceptual  model  of  random 
inhomogeneous  media  as  that  of  a  series  of  anisotropic  Rayleigh  scatterers. 


Because  materials  description  at  mesoscopic  scales  involves  their  specific 
polarimetric  properties,  one  has  to  go  beyond  conventional  effective  medium  approaches. 
For  instance,  one  can  consider  the  mesoscopic  interaction  volume  as  being  a  single 
anisotropic  Rayleigh  scatterer,  defined  by  its  scattering  strength,  degree  of  anisotropy, 
and  orientation  (Fig.  6).  Of  course,  the  properties  of  inhomogeneous  media  are  scale 
dependent  and,  in  general,  the  parameters  of  these  “anisotropic  scatterers”  will  also 
depend  on  the  mesoscopic  volume  considered.  That  means  that  this  local  anisotropy  is 
not  an  intrinsic  property  of  medium;  rather,  it  is  determined  by  the  characteristic  scale 
over  which  the  measurement  is  performed.  The  notion  of  scale-dependent  properties  is  a 


1  For  large  scale  randomly  inhomogeneous  media,  a  multiscale  description  of  the  light  propagation  can  be 
envisioned  where  the  Maxwell’s  equations,  the  transport  equation,  and  the  diffusion  equation  can  be 
applied  to  describe  the  microscopic,  mesoscopic,  and  macroscopic  scale,  respectively. 


powerful  concept,  as  it  provides  novel  means  to  characterize  inhomogeneous  media2. 
This  characteristic  length  scale  of  interaction  may  be  determined,  for  example,  by  the 
particles  dimension  if  we  deal  with  small  inhomogeneous  objects  or  by  the  excitation 
volume,  as  in  the  case  of  near-field  optical  microscopy. 


2.1  Local  anisotropy  of  inhomogeneous  media 

Due  to  their  shape  asymmetries,  small  and  intrinsically  homogeneous  particles  can 
acquire  tensorial  properties  (form  birefringence).  Inhomogeneous  objects  on  the  other 
hand  can  have  optically  anisotropic  properties  even  if  their  components  do  not  have 
shape  asymmetry.  Such  structurally  induced  anisotropy  has  been  observed,  for  instance, 
in  small  spheres  with  eccentric  inclusions.  In  general,  collective  scattering  from 
otherwise  isotropic  scattering  centers  can  manifest  noticeable  anisotropy  in  the  effective 
polarizability. 

Let  us  examine  in  detail  the  mesoscopic  problem  of  an  optically  inhomogeneous 
sphere  with  dimensions  smaller  than  the  wavelength.  The  random  inhomogeneties  inside 
the  sphere  effectively  produce  an  anisotropic  polarizability;  one  method  that  can  be 
applied  in  order  to  measure  this  anisotropic  polarizability  is  that  of  stochastic  scattering 
polarimetry  (SSP).  SSP  will  reveal  the  effective  anisotropic  polarizability  of  an 
inhomogeneous  sphere  if  the  magnitudes  of  the  polarizability  elements  remain  constant 
and  polarimetric  observations  are  performed  for  random  orientations  and  different 
packings. 

It  is  obvious  that  for  some  particular  inhomogeneous  sphere,  the  composite 
particulates  may  be  distributed  resulting  in  an  isotropic  scattering,  while  for  others  of 
similar  properties  it  can  be  very  anisotropic.  To  account  for  situations  where  the 
magnitude  and  orientation  of  the  polarizability  varies,  a  method  can  be  outlined  that  fully 
measures  the  symmetric  polarizability  tensor  of  some  object.  Currently,  in  optics,  one 
can  only  measure  intensities;  therefore  a  complete  determination  of  the  polarizability 
tensor  would  require  three  independent  excitation  fields  and  two  separate  polarimetric 
detectors.  Once  the  complete  symmetric  polarizability  tensor  is  retrieved,  the  diagonal 
form  of  the  polarizability  and  the  corresponding  angles  of  rotation  are  found  using  an 
eigenvalue  decomposition. 

Using  the  above  described  method,  the  diagonal  elements  of  polarizability  tensor 
and  the  angles  specifying  its  orientation  can  be  determined  for  each  realization  of  the 
random  distribution  of  inclusions.  The  situation  can  be  modeled  using  the  coupled  dipole 
approximation  (CDA)  described  in  the  previous  report,  where  dipoles  with  different 
polarizabilities  a  are  randomly  placed  within  the  volume  of  the  sphere.  A  statistically 
relevant  number  of  random  realizations  of  the  sample  should  be  performed.  During  every 
reconstruction  procedure  diagonal  components  of  polarizability  tensor  were  ordered  in 
such  a  way  that  a>b>c.  The  probability  distributions  of  polarizability  components  and 
orientation  angles  are  shown  in  Figure  7.  One  can  see  that  the  distribution  of  angles  is 


2  This  concept  has  been  used  to  describe  different  statistical  regimes  in  near-field  scattering  from  random 
media. 


uniform.  This  behavior  is  specified  by  the  symmetry  of  the  problem:  all  the  orientations 
of  effective  anisotropic  polarizability  of  mesoscopic  sphere  should  be  equally  probable. 
Figure  7  also  shows  that  inhomogeneous  objects  indeed  possess  anisotropic  properties  at 
mesoscopic  scales. 


Figure  7.  Probability  density  functions  of  ordered  diagonal  elements  of 
polarizability  tensor  as  obtained  through  the  use  of  a  sequential  three 
independent  excitations  and  two  polarimetric  detection  schemes. 

Let  us  now  examine  how  the  elements  of  the  polarizability  tensor  can  be  measured.  Of 
course,  a  complete  polarimetric  measurement  involving  three  different  excitations 
followed  by  detection  in  three  different  states  of  polarization  can  provide  all  the 
necessary  information.  However,  in  most  scattering  situations  it  is  difficult  (or  not 
desired)  to  excite  with  three  different  polarizations.  The  question  is,  can  some 
information  can  still  be  retrieved?  In  fact,  if  some  assumptions  about  the  material  can  be 
made,  then  the  anisotropic  polariability  distributions  can  be  found  by  using  one  single 
polarized  excitation,  and  a  single  polarimetrically  measured  intensity. 

Due  to  the  random  nature  of  inhomogeneous  media,  there  is  no  preferential 
structuring  of  the  anisotropic  polarizability  tensor.  That  means  that  the  main  axes  of 
polarizability  tensor  are  oriented  uniformly  in  space.  From  symmetry  considerations  it 
should  also  follow  that  the  values  of  each  of  the  three  diagonal  elements  describing  an 
anisotropic  polarizability  are  random  variables  having  the  same  probability  distributions 
fa  (x)  =  fb  (x)  =  fc  (x) .  If  the  mesoscopic  region  is  excited  by  a  wave  linearly  polarized 
along  x,  then  the  scattered  field’s  magnitude  along  the  same  direction  is: 

Ex  =  vl(0,(p,y/)a  +  v2(O,(p,y/)b  +  v3(0,(p,y/)c, 

Vj  (0,  (p,  Y)  +  v2  (0,  (p,  y/)  +  v3  (0,  (p,  y/)  =  1, 

where  the  functions  v1?  v2  and  v3  describe  the  orientation  of  the  anisotropic 

polarizability  and  are  functions  of  the  random  rotations  0,  cp,  and  y/  .  It  is  possible  to 
show  that  the  random  variables  vt,  v2  and  v3  have  the  same  probability  distribution 
functions  that  can  be  described  by  function 
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However,  these  variables  are  not  independent,  and  their  joint  distribution  may  be  found  to 
be 
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(15) 


Solving  Eq.  (13)  for  unknown  PDFs  fa(x),  fh  (x) ,  /c(v)by  using  measured 

PDF  of  x-polarized  electric  field  and  then  using  the  order  statistics  allows  recovering  the 
PDFs  of  diagonal  polarizability  components.  The  distributions  of  the  ordered  diagonal 
tensor  elements  fuk\(x)  (k  =  1,2,3),  can  be  inferred  from  the  original  probability 

distribution  function  f(x )  and  the  cumulative  distribution  function  F(x) 


fx{k) (x)  =  ^  _  ky  F(Xr  (l -  F(x))n~k  f(x),  n  —  3.  (16) 

For  example,  applying  such  ordered  statistics  constraint  to  a  Gaussian  distribution  with 
mean  /u  and  standard  deviation  a  leads  to  three  other  “Gaussian  like”  distributions  with 
means  //  and  jii±3cj I(2-Jk)  .  This  ordering  allows  us  to  recover  the  three  distributions  of 
the  diagonal  elements  of  the  polarizability,  reducing  the  demand  on  the  measurement: 
single  excitation  polarization  and  single  co-polarized  measurement. 

Our  CDA  simulations  show  that  PDFs  of  a,  b,  and  c  for  random  media  in  many 
cases  may  be  described  by  Gaussian  distribution.  Therefore,  one  can  write: 


fa(a)  =  fb(b)  =  fc(c)  = 
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Resulting  co-polarized  intensity  distributions  are  also  in  many  cases  Gaussian.  Using  the 
assumptions  about  Gaussian  form  of  PDFs,  we  can  solve  numerically  Eq.  (13)  and  find  a 
relation  between  the  moments  of  the  scattered  intensity  distribution  p(\Ex  I2)  and  the 
moments  of  the  Gaussian  distribution  of  polarizability  elements: 
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These  expressions  are  found  from  the  assumptions  that  refractive  index  of  components  of 
inhomogeneous  media  is  bounded  (usually  between  1  and  3).  It  follows  that  the  possible 
variations  of  both  mean  of  intensity  and  standard  deviation  are  also  limited.  Thus,  the 

functions  /  =  f  ( pa ,  cra ),  la  =  fa  ( pa ,  cra )  can  be  decomposed  into  series  with 
respect  to  variables  jua ,  cra.  It  seems  that  a  decomposition  up  to  the  second  order  produce 
very  good  results  as  illustrated  in  Figs  8  and  9  where  two  media  are  considered  that  have 


the  same  effective  properties  (same  effective  refractive  index)  but  different  local 
morphologies  as  can  be  seen.  The  set  of  different  PDF’s  for  polarizibility  elements  are 
different  as  observed  in  Figure  9  indicating  that  a  polarimetric  measurement  can  actually 
discriminate  between  these  two  media. 
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Figure  8.  Example  of  two  random  media  having  the  same  effective  refractive 
index  as  estimated  from  the  Brugememan’s  effective  medium  theory. 
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Figure  9.  Probability  density  functions  of  ordered  diagonal  elements  of  polarizability 
tensors  for  samples  A  and  B  in  Fig.  8  as  obtained  through  the  use  a  single 
polarimetric  excitation. 


2.2  Scale  dependent  local  anisotropic  polarizability 


To  quantitatively  characterize  the  local  anisotropic  polarizability  (LAP),  we 
introduce  an  anisotropy  factor  defined  as  the  contrast  calculated  for  diagonal  components 
of  polarizability  tensor  a 


/3Tr(V)-Tr(tf)2  /-  Jaa (an  -ab)  +  ah ( ah  - a. )  +  ac (a.  - a„ ) 


Tr(«) 


aa  +  ab  +  ac 


where  Tr  denotes  the  trace  of  tensor  a  Note  that,  in  the  past,  other  definitions  have 


been  used  for  such  anisotropy  factor.  In  Ref.  [Error!  Bookmark  not  defined.]  for 
instance,  an  anisotropy  factor  S  was  defined  as  the  variance  of  depolarization  factors 
v(r)  —  \  I  —  f  G0(r,  r'jc/V .  Here  G0(r,r')  is  the  regular  part  of  the  quasistatic 
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free-space  dyadic  Green’s  function  for  electric  field,  /  is  the  unity  tensor.  In  this 
designation,  the  local  anisotropy  factor  cannot  depend  on  the  excitation  volume  and, 
moreover,  its  locality  can  be  violated  in  3D  random  composites  without  structural  self¬ 
similarity,  i.e.  in  nonfractal  composites  that  are  of  interest  here..  The  definition  of  A  in 
Eq.  (19)  is  most  appropriate  for  our  discussion,  which  focuses  on  describing  the  form 
anisotropy  and  not  necessarily  on  the  absolute  magnitude  of  a  specific  dipole  moment. 

Using  the  definition  in  Eq.  (19),  the  anisotropy  factor  A  was  calculated  for  every 
realization  of  the  localized  inhomogeneous  volume.  Of  course,  an  averaged  A  can  then 
be  calculated  from  the  recovered  ensemble  of  values  of  this  parameter.  The  average 
anisotropy  factors  calculated  for  the  two  materials  illustrated  in  Fig.  8  are  9.2- 10~3  and 
13.2-10  3 ,  respectively.  The  30%  difference  clearly  demonstrates  that  LAP  is  a  parameter 
that  can  be  used  to  quantify  differences  between  macroscopically  similar  media. 


In  the  preceding  discussion,  LAP  was  examined  over  one  single  length  scale.  This 
situation  corresponds  to  fixed  volume  of  light-matter  interaction  imposed  by  the 
measurement  procedure.  In  case  of  the  two  different  media  presented  in  Figure  8,  the 
differences  will,  of  course,  diminish  as  this  volume  of  interaction  increases;  the  two 
different  optical  responses  will  converge  toward  the  same  macroscopic  value 
corresponding  to  an  isotropic  polarizability  tensor.  The  rate  of  this  convergence  however 
may  be  different  depending  on  the  specific  structural  morphology. 

We  will  turn  now  our  attention  to  LAP’s  dependence  on  the  volume  of  interaction. 
We  have  repeated  the  previous  analysis  for  spheres  of  different  radii  R  and  the  results  are 
presented  in  Figure  10  where  we  compare  the  anisotropy  factor  A  for  the  case  of  two 
different  sizes  of  spherical  inclusions  randomly  distributed  within  probing  volume  of 
different  sizes.  The  main  observation  is  that  A  always  attains  a  maximum  that  defines  a 
new  length  scale  characterizing  the  electromagnetic  interaction.  This  maximum 
anisotropy  length  (MAL)  represents  the  length  scale  over  which  the  response  of 
inhomogeneous  medium  is  most  sensitive  to  the  polarization  (vector)  properties  of  the 
excitation  field.  In  other  words,  it  is  at  this  scale  that,  in  average,  the  depolarization  of 


light  during  scattering  occurs  more  effectively.  Along  with  scattering  mean  free  path  and 
transport  mean  free  path  that  describe  the  way  the  energy  is  transferred,  MAL  represents 
another  interaction- specific  length  scale  that  characterizes  the  propagation  of  polarized 
fields  through  random  media. 

At  length  scales  smaller  than  MAL,  the  local  polarizability  becomes  isotropic  with 
A  tending  to  zero  as  the  probing  volume  decreases.  At  larger  scales,  the  macroscopic 
behavior  is  gradually  approached  leading  again  to  an  effectively  isotropic  scattering 
volume  with  A  =  0 . 


Figure  10.  Effective  anisotropy  factor  A  as  a  function  of  excitation  volume  R  normalized  by  inclusion 
diameter  d  for  spherical  inclusions  with  a  refractive  index  of  1 .5  randomly  distributed  in  vacuum.  The 
continuous  lines  correspond  to  inclusions  with  diameter  A/ 32  while  the  dashed  lines  correspond  to 

inclusions  of  diameter  d/64 .  Curves  1  to  4  correspond  to  volume  fractions  of  inclusions  of  0.3,  0.2,  0.1, 
and  0.025,  respectively. 


As  can  be  seen  in  Figure  10,  the  values  of  the  anisotropy  factor  appear  to  be 
independent  of  the  size  of  inhomogeneties.  This  happens,  because,  in  our  example,  the 
interaction  inside  the  inhomogeneous  volume  considered  is  mostly  within  the 
electrostatic  regime.  Therefore,  the  behavior  of  A  does  not  depend  on  the  wavelength  and 
is  fully  scalable  with  inclusions’  dimensions.  Also  noticeable  in  Figure  10  is  the  faster 
decay  of  A  for  higher  volume  fractions  of  inclusions  inside  the  sphere  of  interaction. 
This  can  also  be  easily  explained  by  realizing  that,  for  a  given  excitation  volume,  the 
larger  number  of  inclusions  corresponding  to  a  higher  volume  fraction  represents  in  fact  a 
more  isotropic  medium. 

In  the  particular  case  when  the  spherical  inclusions  can  be  considered  as  packed 
hard  spheres,  we  found  that  MAL  has  a  simple  interpretation.  As  illustrated  in  11,  in  this 
case  MAL  defines  the  volume  containing,  in  average,  three  inclusions.  Note  that  three 
inclusions  represent  the  minimum  number  of  particles  necessary  to  form  a  fully 
anisotropic  scatterer.  Thus,  the  statistical  averaging  for  scatterers  containing  more  than 
three  inclusions  results  in  a  gradual  decrease  of  the  anisotropy  factor. 


Figure  11.  Maximum  anisotropy  length  (A)  normalized  by  the  diameter  of  inclusions  d  versus  the  volume 
fraction  of  inclusions  /.  Open  circles  and  crosses  represent  MAL  values  corresponding  to  inhomogeneous 

media  with  inclusions  of  diameters  of  /t/64  ,  zl/32  respectively.  The  solid  line  corresponds  to  the  volume 
containing  on  average  3  inclusions. 

Of  course,  the  other  factors  determining  the  optical  response  of  a  composite  medium  are 
the  intrinsic  properties  of  the  components.  It  is  expected  that,  in  general,  higher 
anisotropy  factors  will  characterize  materials  with  increasing  dielectric  contrasts.  This  is 
evident  in 

Figure  12  were  we  plot  the  value  of  the  maximum  anisotropy  Aniax  as  a  function 
of  dielectric  contrast  of  inclusions  The  calculations  also  demonstrate  that  the  values  of 
Arrax  simply  scale  with  the  magnitudes  of  dielectric  functions  of  components  indicating 

that  the  Aniax  behavior  is  determined  only  by  the  material’s  structure  and  not  by  its 
composition. 


Figure  12.  Maximum  of  anisotropy  factor  Amax  as  a  function  of  volume  fraction/of 
inclusions  with  50nm  in  diameter  and  having  different  refractive  index  contrasts. 


To  conclude,  we  have  demonstrated  that  at  mesoscopic  scales,  the  optical 
response  of  random  media  consisting  of  optically  isotropic  components  may  be 
interpreted  in  terms  of  local  anisotropic  polarizabilities  (LAP).  At  mesoscopic  scales, 
different  materials  can  be  characterized  by  their  specific  anisotropic  polarizabilities  even 
though  they  may  have  similar  effective  dielectric  permittivities  when  described  in  terms 
of  an  effective  medium  approach. 

We  have  also  shown  that  material  properties  at  mesoscopic  scales  depend  on  the 
volume  of  interaction,  in  this  respect,  LAP  is  an  attribute  of  the  electromagnetic  field- 
matter  interaction.  However,  a  characteristic  length  scale,  maximum  anisotropy  length 
(MAL),  exists  at  which  the  degree  of  local  anisotropy  A  reaches  its  maximum.  At  this 
scale  the  inhomogeneous  materials  are  most  sensitive  to  the  polarization  of  incident  light. 
Thus,  electromagnetic  wave  interaction  on  this  scale  length  results  in  the  maximal 
depolarization.  Along  with  other  characteristics  length  scales  such  as  the  scattering  mean 
free  path,  the  value  of  MAL  reflects  essential  intrinsic  properties  of  random  media. 

There  are,  of  course,  different  means  for  characterizing  the  local  optical  properties 
of  inhomogeneous  media.  For  instance,  the  scalar  contrast  of  scattered  intensities 
measures  the  relative  variations  of  the  scattering  cross-section  within  the  interaction 
volume.  When  this  volume  increases,  the  scattered  intensity  variations  decay 
monotonically  to  zero  with  a  rate  depending  on  the  medium’s  properties.  In  this  case 
however,  only  asymptotic  scales  can  be  determined  which  may  affect  the  specificity. 
MAL  on  the  other  hand  is  not  only  derived  from  a  tensorial  feature  of  the  material  but  it 
is  also  a  local  property.  Its  value  is  a  basic  characteristic  of  material’s  morphology. 

Our  findings  are  relevant  to  the  design  of  novel  materials  because  this  new 
electromagnetic  interaction  scale  represents  the  material  scale  at  which  the  polarimetric 
response  of  a  medium  is  most  sensitive  to  the  excitation  field. 


3.  Field  transfer  matrix  of  scattering  media 


Multiple  scattering  is  typically  considered  to  degrade  the  information  in  a  beam 
propagating  through  a  random  medium.  For  temporally  coherent  radiation  propagating 
through  an  optically  thick  medium,  the  light  scattered  out  of  the  medium  will  produce  an 
interference  pattern  with  alternating  bright  and  dark  regions  known  as  speckles.  A 
speckle  pattern  having  a  negative  exponential  intensity  distribution  is  called  fully 
developed  and  obeys  Gaussian  field  statistics.  Since  the  field  distribution  is  a  Gaussian 
random  variable,  it  was  thought  to  contain  neither  information  about  the  underlying 
scattering  medium  nor  information  about  the  beam  incident  on  the  medium. 

We  investigated  a  new  approach  to  recovering  information  from  a  speckle  pattern 
that  does  not  rely  on  ensemble  averaging  or  assumptions  about  the  underlying  field 
statistics.  We  have  actually  measured  the  components  of  the  field  transfer  matrix  across 
a  particular  output  plane.  In  a  proof  of  concept  demonstration,  we  show  successful 
recovery  of  the  polarization  state  of  an  unknown  beam  from  the  speckle  pattern  it 
produces.  Our  technique  is  a  deterministic  procedure  performed  simultaneously  in  many 
different  spatial  locations.  Because  the  technique  does  not  rely  on  any  assumptions  about 
the  statistics  of  the  field  distribution  and  because  it  is  carried  out  by  point  operations 
rather  than  image  correlations,  it  can  be  used  in  regimes  ranging  from  no  scattering  to 
high-order  multiple  scattering. 


Consider  a  linear,  multiply  scattering  medium  illuminated  by  polarized  light.  The 
field  will  undergo  many  series  of  scattering  events  before  emerging  from  the  medium.  In 
a  particular  series,  the  field  after  the  first  scattering  event  at  iq  that  illuminates  the  second 
scatterer  at  position  r2  may  be  expressed  as 

E(r2  =  «(r2,r,  )E(ri  )inc  (20) 

where  a  is  a  complex  tensor  that  expresses  the  magnitude  and  phase  of  the  coupling 
between  incident  and  scattered  field  components.  It  also  includes  transferring  the 
scattered  field  to  the  point  r2.  Now,  replacing  Einc  with  Escat,i  and  utr,,!-, )  with  a(r3,r2) , 
we  obtain  the  field  at  position  r3  on  the  third  scatterer.  This  process  can  be  continued 
until  the  field  reaches  the  detector  at  position  r  to  obtain  the  contribution  of  a  particular 
scattering  path  to  the  field  at  the  detector: 


E(r)det  = 


n«(ri+i>ri) 


E(iq) 

inc  (t(r,  iq  )ejJ  Efiq  )inc 


(21) 


In  this  particular  decomposition  the  scattering  matrix  that  relates  the  input  and  output 
field  contains  information  about  the  propagation  from  rj  to  iq+i  as  well  as  information 
about  the  scattering  event  at  rj.  Since  we  are  interested  in  the  transfer  of  the  field  through 
the  medium,  not  in  the  particular  details  of  how  the  transfer  occurs,  the  mixing  of 
scattering  and  propagation  information  is  not  important. 

Furthermore,  there  are  many  scattering  paths  connecting  the  illuminated  points  on  the 
front  surface  of  a  multiply  scattering  medium  to  each  detection  point  after  the  medium  as 
shown  in  Fig.  13.  If  the  scattering  path  lengths  are  shorter  than  the  coherence  length  of 
the  illuminating  radiation,  the  scattering  process  will  not  depolarize  the  output  light;  it 


will  only  change  the  radiation’s  state  of  polarization.  The  resulting  field  at  any  detection 
point  due  to  a  particular  input  point  is  simply  the  coherent  summation  of  the  outputs  of 
each  of  the  different  paths  starting  at  the  illuminated  point  on  the  front  surface  and  ending 
at  the  detection  point.  Also  shown  in  Fig.  13  are  three  random  walks  on  the  Poincare 
sphere  that  begin  in  the  same  polarization  state  (chosen  here  to  be  circular).  However, 
since  they  interact  with  different  configurations  of  the  scattering  medium,  they  arrive  at 
the  detector  in  different  states  of  polarization,  denoted  by  the  large  dots  of  the  appropriate 
color.  The  total  field  at  the  detector  is  the  coherent  sum  of  the  different  paths  and  is  in 
the  state  of  polarization  marked  by  a  white  dot  labeled  “detected”  on  the  Poincare  sphere. 
Thus,  the  resulting  field,  neglecting  time  dependence,  at  a  point  r  due  to  all  illuminated 
input  points  can  be  written  as 

E(r ) total  =  X,  Yjeff  “(r’  r'>  W  E(rn  )fac  =  X„  “(r’  F"  )E(r»  c  =  (22) 

The  normalized  input  field,  e,-nc,  can  be  factored  out  of  the  sum  so  long  as  its 
polarization  is  spatially  constant,  and  we  can  express  the  resulting  output  field  in  terms  of 
a  single  transfer  matrix,  a(r) ,  which  subsumes  the  intensity  profile  of  the  illumination. 
For  a  given  experimental  geometry  and  illumination  source,  a(r)  is  only  a  function  of 
detector  location;  however,  if  the  spatial  intensity  profile  of  the  illumination  is  modified, 
«(r)  will  change  as  well  because  the  intensity  profile  acts  as  a  weighting  function  for  the 
contribution  of  each  scattering  path.  The  tilded  quantities  represent  the  only  measurable 
parameters  of  the  scattering  material  because  it  is  not  practically  possible  to  separate  the 
contributions  of  individual  paths  from  the  detected  intensity.  It  is  also  important  to  note 
that  the  illuminated  points  can  have  an  arbitrary  spatial  extent  and  spatial  intensity  profile 
on  the  random  scattering  medium  as  long  as  the  scattering  paths  remain  coherent  with 
one  another. 


Fig.  13.  Random  walks  through  a  static  random  medium  and  their  resulting 
change  in  polarization  state. 


The  problem  of  determining  the  elements  of  the  transfer  matrix  to  a  given  point  in  the 
detection  plane  does  not  depend  on  the  precise  nature  of  the  process  that  produced  the 
output  field  and  thus  applies  to  all  scattering  regimes.  In  general,  a  has  nine  elements 
with  eighteen  unknowns:  nine  coupling  magnitudes  and  nine  phases.  Because  the 
scattering  is  not  isotropic,  a  depends  on  the  direction  of  propagation  of  the  incident  light, 
and  it  is  difficult  to  measure  the  full  transfer  matrix  for  an  arbitrary  geometry.  However, 


if  the  scattering  medium  is  surrounded  by  an  isotropic  medium,  the  electric  field  of  the 
illumination  is  confined  to  a  plane  and  can  be  decomposed  into  two  orthogonal 
polarization  states  with  a  phase  between  them.  If  the  scattered  fields  are  allowed  to 
propagate  away  from  the  scattering  medium  before  detection,  they  can  also  be 
decomposed  into  two  orthogonal  polarization  states  and  a  phase  term.  In  this  situation, 
the  transfer  matrix  has  only  four  elements  consisting  of  eight  unknowns,  which  can  be 
determined  by  illuminating  the  scattering  medium  with  appropriately  polarized  light. 
Moreover,  we  can  choose  one  of  the  elements  of  the  transfer  matrix  to  be  real  since  we 
cannot  measure  absolute  phase  at  a  point  and  only  compare  the  intensities  between 
points.  Additional  simplifications  can  be  introduced  by  realizing  that  it  is  not  necessary 
to  characterize  all  seven  of  the  remaining  unknowns  simultaneously.  If  a  polarizer 
oriented  along  the  x-axis  is  placed  between  the  scattering  medium  and  the  detector, 
a2l=a22=0,  and  there  are  only  three  unknowns  that  need  to  be  characterized.  The 
detected  intensity  at  point  r  is  then  given  by 

I(r)  =  [au(r)Exf  +  [a{2(r)Eyf  +  2au(v)an(v)ExEy  cos[6>  +  ^(r)] ,  (23) 

where  0  is  the  phase  between  the  x  and  y  components  of  the  incident  field,  and  Ex  and  Ey 
are  their  respective  magnitudes.  In  Eq.  (23),  (p  is  the  phase  introduced  by  the  coupling  of 
Ey  into  a  scattered  x  polarized  field,  and  «, ,  and  «12  are  the  magnitudes  of  the  coupling  of 
the  incident  x  and  y  polarized  fields,  respectively,  into  scattered  fields  polarized  along  x. 
The  polarization  of  the  scattered  field  and  the  elements  of  fi  measured  are  determined  by 
the  orientation  of  the  final  polarizer. 


3.1  Calibration  procedure 

The  procedure  for  recovering  an  unknown  incident  state  of  polarization  comprises  three 
main  steps: 

•  measuring  the  transfer  matrix  for  many  points  in  the  detection  plane, 

•  selecting  unique  field  combinations  or  transfer  matrices,  and 

•  solving  Eq.  (10)  for  the  incident  field  components  using  the  transfer  matrices. 

First,  the  needed  components  of  the  transfer  matrix  corresponding  to  each  detector  point 
are  determined  via  a  calibration  with  known  fields.  To  determine  the  magnitudes  of  the 
elements  of  the  transfer  matrix  and  eliminate  the  sign  ambiguity  in  the  argument  of  the 
cosine  in  Eq.  (10),  we  use  four  calibration  states.  The  relative  spatial  intensity  profiles  of 
the  unknown  source  and  the  calibration  source  should  be  the  same  since  the  intensity 
profile  of  the  illumination  weights  the  contribution  of  each  scattering  path  to  the  detected 
intensity. 

A  fully  polarized  field  is  characterized  by  three  different  parameters,  and,  as  with  the 
transfer  matrix  elements,  at  least  three  independent  combinations  of  the  incident  field  are 
needed  to  completely  determine  its  polarization  state.  For  a  static  system,  each  point  in 
the  detection  plane  sees  a  particular  transformation  of  the  input  field  resulting  from  the 
combination  of  the  scattering  paths  that  end  at  that  point.  Thus,  determining  the 
parameters  of  a  fully  polarized  field  requires  that  the  intensity  be  measured  at  three  or 


more  points  with  independent  transfer  matrices.  As  a  result,  the  detection  system  must 
resolve  at  least  three  speckles. 

As  an  example,  the  use  of  a  particular  combination  of  detectors  is  illustrated  in  Fig. 
14.  The  axes  of  this  representation  of  transfer  matrices  are  defined  in  a  manner 
analogous  to  the  Poincare  sphere  with  the  elements  of  the  transfer  matrices  at  a  particular 
point  taking  the  place  of  the  input  field  components  that  they  couple 

A  =  (an2  -5122)/(5n2  +  5122) 

cr  =  2al !«!  2  cos(^ )  l(d{  2  +  ax  2 )  (24) 

(f)  =  25,  ,5, 2  sin(<^ )  /(5,  2  +  ax  2 ) 

The  points  shown  in  blue  have  a  transfer  matrix  of  an  *  0  and  dl2~  0 ,  while  the  points 
shown  in  green  have  a  transfer  matrix  of  5, , «  0  and  dl2  *0 .  Because  the  output  points 
represented  by  the  blue  and  green  areas  couple  only  one  of  the  two  input  field 
components,  they  measure  the  x  and  y  components  of  the  unknown  field.  The  points 
shown  in  red  on  the  other  hand  have  au  «  a12  and  contain  the  information  about  the  phase 
of  the  unknown  field  because  both  of  the  input  field  components  are  coupled  into  the 
measured  intensity. 


Fig.  14.  Example  of  groups  of  transfer  matrices  that  can  be  used  to  recover  an 
unknown  incident  field. 


Since  there  is  no  reason  to  select  any  particular  group  of  independent  detection  points, 
it  is  possible  to  form  many  different  groups  and  then  perform  a  statistical  analysis  on  the 
recovered  fields  rather  than  relying  on  the  result  of  a  single  combination.  Note  that  a 
medium  with  transfer  matrices  covering  only  a  small  portion  of  the  sphere  may  not  be 
used  to  fully  analyze  any  unknown  field  from  its  speckle  pattern  simply  because  the 
medium  does  not  produce  a  sufficient  number  of  independent  combinations  of  the  field. 

Finally,  for  the  selected  groups  of  transfer  matrices  as  illustrated  in  Fig.  2,  the 
unknown  field  parameters,  Ex,  Ey,  and  9 ,  are  determined  by  solving  the  following  system 
of  equations  for  each  group  of  detectors 

7(r1)  =  [511(r1)£x]2 
/(r2)  =  [512(r2)^]2 

I(r3)  =  [an(r3)Ex]2  +[512(r3)E>,]2  +  2511(r3)512(r3)E*Eycos[<9  +  ^(r3)] 


(25) 


where  all  of  the  transfer  matrix  elements  are  known  from  the  calibration  process  and 
where  rj  denotes  the  location  of  the  point  in  the  detection  plane  rather  than  the  location  of 
a  scatterer. 

3.2  Proof  of  concept  demonstration 

In  order  to  demonstrate  that  our  procedure  is  valid  across  all  scattering  regimes,  we 
measured  transfer  matrices  for  the  extreme  cases  of  deterministic,  single  scatterers  and  a 
heavy  multiple  scatterer.  For  the  single  scatterers,  the  measured  transfer  matrices  can  be 
qualitatively  compared  to  the  expected  values  using  the  spherical  representation  from  Fig. 
14.  On  the  other  hand,  since  the  multiply  scattering  medium  produces  many  diverse 
mixings  of  the  incident  field,  it  is  difficult  to  assess  the  capability  of  our  procedure  by 
direct  examination  of  the  measured  transfer  matrices.  In  this  case,  we  will  use  the 
transfer  matrices  to  recover  the  states  of  polarization  of  plane  waves  in  order  to  show  our 
procedure  works  even  in  a  regime  of  heavy  multiple  scattering. 

In  our  experiment,  a  laser  beam  with  a  controlled  state  of  polarization  was  incident  on 
the  scattering  medium  and  the  speckle  pattern  resulting  from  its  transmission  through  the 
medium  was  recorded  by  a  CCD  camera.  A  polarizer  with  a  fixed  orientation  was  placed 
in  front  of  the  detector  to  simplify  the  analysis  as  explained  before.  The  speckles 
produced  by  the  multiply  scattering  medium  were  approximately  five  pixels  across  on  the 
CCD  and  each  measurement  sampled  a  few  thousand  speckles. 

Figure  15  illustrates  the  measured  transfer  matrices  for  a  polarizer,  a  quarter- wave 
plate,  and  a  multiply  scattering  medium.  The  multiply  scattering  medium  is  a  composite 
dielectric  material  with  a  thickness  of  100pm  and  characterized  by  a  transport  mean  free 
path  of  10pm.  Figures  (a)  and  (b)  show  measurements  from  a  polarizer  oriented  at 
approximately  ±45°  and  a  quarter-wave  plate  rotated  in  15°  increments  from  0°  to  90°. 
From  linear  optics  theory,  we  expect  a  polarizer  oriented  at  45°  and  135°  to  have  equal 
coupling  strengths  through  the  final  polarizer  for  both  incident  x  and  y  field  components. 
Also,  when  the  polarizer  is  oriented  at  45°,  the  transmitted  field  components  will  be  in 
phase,  and  when  it  is  oriented  at  135°,  the  transmitted  field  components  should  have  a  n 
phase  difference  between  them.  In  (a)  we  see  comparable  coupling  of  the  orthogonal 
field  components,  although  there  is  a  slight  misalignment  of  the  polarizer,  and  rotating 
the  polarizer  from  45°  to  135°  introduces  a  ji  phase  shift  between  the  elements  of  the 
transfer  matrix  as  expected.  In  (b),  the  green  line  denotes  the  path  that  the  transfer  matrix 
of  a  quarter-wave  plate  should  follow  on  the  sphere  as  the  wave  plate  is  rotated  through 
90°,  and  the  labeled  groups  of  points  are  the  measured  transfer  matrices  for  a  quarter- 
wave  plate  in  the  indicated  orientations.  In  both  cases  the  measured  transfer  matrix 
elements  show  behavior  consistent  with  what  was  expected  and  demonstrate  the  ability  of 
our  process  to  analyze  deterministic  single  scatterers. 


Fig.  15.  Effective  transfer  matrices  measured  for  (a)  A  polarizer  oriented  at 
roughly  45°  and  135°.  (b)  A  quarter-wave  plate  rotated  by  90°  in  15°  increments. 

(c)  A  multiply  scattering  solid  sample. 

In  Fig.  15(c)  we  present  a  subset  of  the  measured  transfer  matrices  for  the  multiply 
scattering  sample.  For  a  truly  random  medium,  we  would  expect  the  measured  transfer 
matrices  to  uniformly  cover  the  sphere  of  possible  matrices;  however,  because  of 
numerical  instabilities  involving  calculations  with  small  numbers,  the  points  in  a  ring 
around  A  axis  of  the  sphere  are  rejected  by  our  processing  algorithm.  Near  the  axis, 
either  an  »  512  or  «  «l2,  and  the  smaller  of  the  two  can  be  approximated  as  0  so  that 
transfer  matrices  very  close  to  the  A  axis  are  moved  onto  the  axis.  Also,  there  seems  to 
be  some  clustering  of  the  points  near  </>=-!.  Even  though  the  transfer  matrices  are  not 
uniformly  distributed  on  the  sphere,  our  measurements  show  that  the  sample  produces  a 
sufficiently  large  number  of  substantially  different  mixings  of  the  incident  field.  This  is 
the  only  requirement  for  recovering  the  state  of  polarization  of  the  incident  field. 

Since  the  medium  exemplified  in  Fig.  15(c)  is  in  a  regime  of  multiple  scattering,  it  is 
difficult  to  assess  the  accuracy  of  the  transfer  matrix  measurement  by  viewing  of  the 
matrices  using  this  spherical  representation.  However,  one  can  still  examine  their 
accuracy  by  using  the  determined  matrices  to  infer  the  polarization  states  of  different 
beams  illuminating  the  medium.  Fig.  16  shows  typical  experimental  results  for  a  +45 
degree  linear  polarization  and  an  elliptical  input  state,  represented  by  the  blue  dots.  The 
white  dots  denote  the  polarization  states  recovered  by  different  pixel  groups,  and  the  red 
dots  represent  the  geometric  centers  of  the  white  data  points. 
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Fig.  16.  Poincare  sphere  representations  of  the  polarization  states  recovered  by 
individual  pixel  groups  (white  dots),  the  geometric  centers  of  the  clouds  of  white 
points  projected  onto  the  surfaces  of  the  spheres  (red  dots),  and  the  expected 
polarization  states  (blue  dots). 

The  experimental  data  shown  above  subtend  a  solid  angle  of  approximately  0.158 
steradians.  The  spread  in  the  states  recovered  by  the  different  detector  combinations  can 
be  mitigated  to  a  large  extent  by  averaging  the  results  of  many  combinations.  The 
averaging  can  be  done  well  because  transfer  matrices  are  measured  for  a  large  number  of 
detectors.  In  order  to  quantify  the  error  in  the  recovered  polarization  state,  we  can 
compare  the  normalized  Stokes  vector  of  the  recovered  field  to  that  of  the  incident  field. 
In  (a)  we  input  a  Stokes  vector  of  (0,  1,  0)  and  measured  (0.018,  0.9997,  0.014),  and  in 
(b)  we  input  (0,  0.643,  0.766)  and  measured  (0.067,  0.638,  0.764).  As  can  be  seen,  the 
recovered  Stokes  vector  components  do  not  deviate  by  more  than  1%  from  their  expected 
value. 

We  have  also  simulated  our  experiment  to  study  the  effects  of  the  detector  selection 
criteria  and  detector  noise  on  the  measured  field  transfer  matrices.  We  choose  5n  and  «l2 
randomly  distributed  uniformly  between  0  and  1  and  to  be  randomly  distributed 
uniformly  between  0  and  2n.  The  magnitudes  of  the  coupling  matrices  were  then  scaled 
so  that  the  resulting  calculated  intensity  distribution  was  similar  to  the  actual  data  for  unit 
strength  electric  field  inputs.  Gaussian  white  noise  with  a  signal  to  noise  ratio  of  34  was 
then  added  to  the  intensity  image  to  simulate  the  detector  noise  in  a  real  measurement. 
Speckle  images  were  generated  for  both  the  calibration  and  test  states  and  processed 
using  the  same  code  as  the  experimental  data.  Our  simulations  indicate  that  the  most 
significant  source  of  error  in  the  data  collection  and  processing  is  the  noise  in  the  detector 
itself.  For  simulated  data  with  no  noise,  the  solid  angle  covered  by  the  measurement  data 
on  the  Poincare  sphere  (white  dots  in  Fig.  4)  is  approximately  0.013  steradians.  When 
noise  comparable  to  the  noise  of  the  detector  used  in  the  experiment  is  added  to  the 
simulated  data  before  processing,  the  spread  of  the  recovered  states  increases  to  0.048 
steradians.  The  remaining  error  is  likely  due  to  mechanical  instabilities  in  the  experiment. 

In  summary,  we  have  demonstrated  that  the  field  transfer  matrices  of  a  system  can  be 
practically  measured  for  specific  geometries.  In  doing  so,  no  assumptions  about  the 
statistics  of  the  random  field  and  no  specific  description  of  the  scattering  process  are 


needed.  We  have  demonstrated  that  a  random  medium  can  be  used  as  a  massively  parallel 
polarimeter  having  a  remarkable  precision.  This  kind  of  polarimeter  should  be  useful  for 
applications  (i)  requiring  “instantaneous”  (one  shot)  or  very  high  speed  measurements, 

(ii)  involving  high  energy  pulses  so  long  as  nonlinear  effects  are  minimal  (this  can  be 
further  investigated),  (iii)  requiring  transfer  of  information  encoded  in  the  polarization 
state  of  light  (i.e.  keying  with  a  large  number  of  available  symbols  per  time  slot). 

Unique  characteristics  of  such  scattering  based  device  are  (i)  wide  acceptance  angle  and 
usable  waveband,  (ii)  the  optical  train  does  not  need  to  be  aligned  with  incident  light,  and 

(iii)  the  massively  parallel  measurement  allows  different  approaches  for  optimizing  the 
analysis. 


4.  Characterizing  the  fluctuations  of  scattered  waves 

The  interaction  of  coherent  electromagnetic  beams  with  most  common  materials  leads  to 
a  variety  of  specific  fluctuations.  As  a  result,  the  electromagnetic  beams  acquire  usually  a 
certain  degree  of  randomness.  The  statistics  of  random  complex  fields  give  rise  to  a 
number  of  measurable  distributions,  including  intensity,  phase  difference,  ellipticity,  and 
states  of  polarization.  The  most  familiar  occurrence  of  random  complex  fields  is  the 
family  of  speckle  phenomena  resulting  from  the  coherent  illumination  of  diffusive 
random  media  and  scattering  from  rough  surfaces. 

A  light  scattering  experiment  in  a  random  medium  provides  a  wealth  of  information 
about  the  properties  of  the  material  under  investigation.  Unfortunately,  the  information  is 
often  too  complex  to  process  in  a  practical  matter.  For  this  reason,  an  ensemble  average 
over  several  material  realizations  is  made  in  order  to  determine  the  mean  statistical 
properties  of  the  sample.  For  example,  coherent  light  that  is  injected  into  a  random 
medium  produces  a  speckle  pattern  that  is  indicative  of  the  medium’s  specific  structural 
properties.  In  order  to  learn  anything  about  the  global  material  properties,  numerous 
speckle  patterns  resulting  from  different  realizations  of  the  material  must  be  averaged. 
This  averaging  inherently  discards  “excess”  information  specific  to  particular 
realizations. 

A  random  medium  is  typically  characterized  by  an  ensemble  of  realizations  of  disorder. 
When  waves  interact  with  a  random  medium,  each  member  of  this  ensemble,  i.e.  each 
particular  realization  of  disorder,  has  its  own  pattern  of  fluctuations  in  the  scattered  wave. 
The  interaction  is  a  non-self-averaging  process,  and  the  complicated  features  of  the 
scattered  waves  are  all  rooted  in  the  structural  properties  of  the  specific  realization  of 
randomness.  In  principle,  the  inverse  problem  can  be  solved  if  (i)  the  phase  coherence  is 
maintained  over  the  entire  interaction,  (ii)  the  process  is  not  dissipative,  and  (iii)  the 
disorder  does  not  vary  in  time.  In  practice,  however,  due  to  finite  sizes  and  experimental 
noise,  one  always  infringes  at  least  on  the  second  requirement.  Furthermore,  the 
information  available  is  often  too  complex  to  process  in  a  practical  manner;  therefore,  an 
average  over  an  ensemble  of  disorder  realizations  is  usually  taken  to  determine  mean 
statistical  properties.  For  instance,  numerous  speckle  patterns  resulting  from  different 
interactions  must  be  averaged  to  learn  about  the  global  properties.  Unfortunately,  this 
averaging  inherently  discards  information  specific  to  particular  realizations  as  well  as  the 
variations  from  one  realization  to  the  next. 

As  such,  the  following  question  arises:  can  one  learn  anything  about  the  stochastic 
process  by  examining  a  set  of  its  individual  realizations?  In  principle,  having  access  to  a 
number  of  samples  of  the  random  process  should  allow  us  to  study  the  rate  of 
convergence  toward  the  ensemble  average  characteristics.  Using  the  example  of 
stochastic  interaction  between  coherent  waves  and  random  media,  we  will  demonstrate 
that  the  significant  sample  to  sample  fluctuations  can  be  used  to  extract  information  not 
available  in  the  ensemble  average.  Because  the  convergence  of  the  statistical  properties 
of  moments  is  a  general  problem  for  numerous  physical  phenomena  that  are  described  as 
random  processes,  concepts  similar  to  the  one  discussed  here  may  be  exploited  in  other 
situations. 


Let  us  consider  the  wave  interaction  with  a  random  medium  characterized  by  a  number 
density  Nv  of  scattering  centers  and  by  the  scattering  cross-section  cr  describing  the 
properties  of  a  single  scattering  event.  For  each  realization  of  disorder  a  ,  the  interaction 
will  be  defined  by  a  specific  distribution  pa  (5)  of  available  wave  path-lengths  s  through 

the  medium.  When  an  ensemble  average  is  taken  over  many  such  realizations,  the  wave 
interaction  will  be  described  by  a  probability  distribution  function 
p(s)=  (pjs))  =  f(s,D )  that  has  a  universal  behavior  depending  only  on  the  normalized 

diffusion  coefficient  D  =  \l(Nvcr) .  Note  that  all  the  experimentally  observable  properties 
of  the  stochastic  interaction  can  be  described  in  terms  of  the  probability  distribution  p(s ) 
whose  exact  functional  form  f(s,  D )  may  also  depend  on  the  particular  geometry  of  an 
experiment.  Clearly,  there  could  be  many  dissimilar  media  with  different  Nv  and  cr  that 

nevertheless  display  the  same  characteristics  upon  ensemble  averaging.  In  practice,  this 
ensemble  can  be  acquired  in  different  ways  for  dynamic  or  stationary  systems,  but  the 
final  result  is  the  same:  the  number  density  and  the  scattering  cross-section  are  being 
coupled  through  the  diffusion  coefficient,  and  only  their  product  is  accessible. 

For  a  specific  realization  of  the  material  disorder,  however,  the  distribution  of  available 
path-lengths  pa(s )  will  deviate  from  the  one  corresponding  to  the  ensemble  average: 

Pa  (S)  =  f(S >  ^I1  +  5 a  (*>  4  )}  (26) 

Because  this  deviation  Sa(s,^a)  is  specific  to  a  particular  realization  of  disorder,  it 
depends  on  variables  not  present  in  the  ensemble  average.  Specifically,  this  can  be 
expressed  through  a  configuration  function  E,a  describing  the  particular  morphology  of 

the  given  realization  a .  The  function  describes  the  locations  of  scattering  centers 

available  in  the  realization  and,  therefore,  depends  only  on  the  number  density  and  not  on 
the  scattering  cross-section.  This  is  the  reason  why  by  examining  the  statistical  properties 
of  8 a  (.v,  8a ),  one  could  infer  information  not  available  in  the  ensemble  average. 

Two  observations  about  the  general  behavior  of  8 a  (.v,  8a )  are  worth  making.  First,  as  the 
length  s  of  the  wave's  path  increases,  more  and  more  different  trajectories  of  the  same 
length  are  possible  through  the  medium,  and  pa(s )  approaches  the  value  corresponding 
to  the  weight  of  trajectories  of  length  s  in  the  ensemble  average.  In  other  words,  in  terms 
of  the  variable  5,  the  function  8a(s,^a)  represents  a  nonstationary  random  process. 
Second,  because  upon  ensemble  averaging  a  scattering  center  will  exist  at  any  position, 
this  random  function  is  of  zero  mean,  (8a(s,<^a))  =  0  (the  mean  is  of  little  practical 

interest).  However,  because  of  the  implicit  dependence  on  Nv ,  it  is  expected  that  higher 
order  statistics  of  Sa(s,ga)  can  be  used  to  reveal  characteristics  of  the  wave-matter 
interaction  not  included  in  the  value  of  D. 


Fig  17.  Sketch  of  path-length  distributions  for  two  media  with  identical  mean 
properties  (same  D).  The  two  media  consist  of  scatterers  of  different  cross- 
sections  a  and  different  number  densities  Nv  and  are  examined  over  the  same 
range  of  path-lengths  .s .  The  medium  with  smaller  number  density  provides 
fewer  possible  paths  of  given  length  s  resulting  in  larger  fluctuations  of  pa  ( s ) . 


This  concept  is  illustrated  schematically  in  Fig.  17  where  the  path-length  distributions 
corresponding  to  two  different  media  are  sketched  over  similar  ranges  of  s .  The  two 
media  consist  of  scattering  centers  having  different  cross-sections  but  also  different 
number  densities  such  that,  upon  ensemble  average,  they  are  described  by  the  same 
diffusion  coefficient.  Clearly,  when  compared  to  all  potential  trajectories,  there  are  fewer 
available  paths  of  given  length  s  through  the  medium  with  less  scattering  centers. 
Consequently,  the  path-length  distribution  pa(s )  deviates  more  significantly  from  the 

ensemble  average  p(s)=  ( pa(s )) .  A  measure  of  the  nonstationary  fluctuations  in  pa{s) 
should  discriminate  between  the  two  media,  as  we  will  show  in  the  following. 

Evidently,  the  random  function  pa  (5)  displays  not  only  fluctuations  in  .v  but  also 


differences  from  one  material  realization  a  to  another.  There  are  many  ways  in  which 
the  two-dimensional  statistical  characteristics  of  pa\s)  can  be  quantified.  Of  course,  a 

simple  averaging  over  a  will  provide  a  path-length  distribution  p(s)=  f(s,D)  which 
corresponds  to  the  ensemble  average.  For  a  single  realization  a  on  the  other  hand,  higher 
order  moments  of  the  fluctuations  in  pa{s )  can  be  evaluated.  Even  though  pa(s)  is 


nonstationary  in  s ,  one  can  still  calculate  simple  estimators  such  as,  for  instance,  the 

1  2 

variance  Va(^a)=  jS^(s,<^a)ds - nSa(s,^a)ds  of  the  fluctuations  along  s.  This  simple 


estimate  however  is  inadequate  because  8a (s, ga)  is  a  zero-mean  random  function  and, 
consequently,  a  unique  and  meaningful  normalization  is  difficult  to  define. 


As  the  deviation  Sa(s,^a)  from  the  ensemble  average  can  be  regarded  as  a  form  of 

disorder,  we  can  choose  to  examine  its  variance  in  terms  of  the  Shannon  information 
entropy: 
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In  Eq.  (27),  we  define  this  finite  scale  entropy  to  account  for  realistic  situations  of  any 
measurement  that  extends  over  a  finite  range  [v, ,  ,s2  ] .  Furthermore,  the  finite  scale 
entropy  can  be  normalized  to  its  maximum  allowable  value  for  the  entire  range 
S  =  s2  -  as 


ha(SVS2)  =  -Ha(SVS2)n°g 


r  1  x 


(28) 


Of  course,  the  normalized  entropy  ha(s1,s2)  will  still  vary  from  realization  to  realization 


and  one  can  further  build  its  average  ha  (Lv, ,  s2 )  over  the  number  of  realizations  available. 
Being  constructed  in  terms  of  the  specific  fluctuations  of  each  realization  a  ,  this  average 
is  a  comprehensive  measure  of  the  overall  fluctuations  in  Sa(s,£a).  It  depends  directly 

on  Nv  through  its  influence  on  the  configuration  function  E,a . 


We  conducted  an  experiment  to  demonstrate  that  structural  information  not  present  in  the 
ensemble  average  can  be  determined  by  analyzing  the  fluctuations  in  realizations  of 
pa(s).  The  distribution  of  photon  path-lengths  through  different  multiply-scattering 

media  was  measured  interferometrically  using  the  procedure  of  optical  path-length 
spectroscopy  (OPS)3.  In  OPS,  the  short  coherence  length  of  radiation  allows  selecting  a 
specific  length  s  of  the  optical  paths  in  a  manner  similar  to  time-of-flight  experiments. 
The  OPS  measurements  are  based  on  fiber  optic  arrangements  that  permit  different 
modalities  for  injecting  light  into  and  collecting  light  from  a  scattering  medium.  The 
configuration  can  be  monostatic,  where  the  same  fiber  acts  as  both  the  source  and 
detector,  or  bistatic  when  the  injection  and  detection  points  are  separated  by  an  adjustable 
distance  A  allowing  for  an  experimental  control  over  the  volume  of  interaction.  In  the 
frame  of  diffusion  theory  for  lossless  media,  the  path-length  distribution  p(s.  A) 
corresponding  to  the  ensemble  average  can  be  evaluated  to  be 
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where  z(,  is  the  so-called  extrapolation  length. 

We  examined  two  different  highly  diffusive  media  that  have  approximately  identical 
average  properties.  These  random  dielectrics  are  essentially  non-absorbing  polymer 
structures  and  have  average  pore  sizes  of  0.45  pm  and  1.2  pm.  Upon  ensemble 
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averaging,  both  are  characterized  by  the  same  value  of  the  transport  mean-free  path  of 
10  ju  m. 
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Figure  19  The  averaged  backscattered  intensities  for  medium  A  (blue  solid 
line)  and  medium  B  (red  dashed  line).  The  insets  show  typical  micrographs  of 
the  materials  examined. 

Path-length  distributions  averaged  over  ten  different  realizations  of  these  random  media 
are  shown  in  Fig.  19  together  with  their  corresponding  scanning  electron  micrographs. 
Even  though  their  structural  morphologies  are  obviously  different,  the  similar  behavior  of 
p(s )  is  a  clear  indication  that,  in  average,  the  two  media  are  being  described  by  the  same 
diffusion  coefficient.  On  the  other  hand,  the  fluctuations  from  the  average  are  rather 
dissimilar  as  can  be  seen  Fig.  20  where  we  plot  the  typical  mean  square  of  the 
fluctuations  8 ?(s,ga)  corresponding  to  the  two  media.  In  general,  the  medium  A  exhibits 
smaller  deviations  from  the  average  which  can  be  interpreted  as  a  larger  number  of 
scattering  trajectories  available  for  each  s .  Note  also  that  the  fluctuations  in  82a  (,v,  ) 

decrease  for  larger  values  of  s  because  these  random  processes  are  nonstationary  as 
discussed  above. 
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Figure  20  Typical  mean  square  fluctuations  rij  (.v.  £,tl )  of  path-length 
distributions  for  media  A  and  B  shown  in  Fig.2. 


Figure  21  Typical  measurment  of  pathlength  distribution  in  a  bistatic 
geometry.  By  increasing  the  source-detector  separation  A  one  effectively 
enlarges  the  volume  of  interaction  with  the  random  medium. 

Let  us  now  examine  in  more  detail  the  situation  where  the  scale  of  available  path-lengths 
is  varied.  In  practice,  this  amounts  to  controlling  the  size  of  the  interaction  volume  which 
can  be  implemented  in  the  bistatic  OPS  measurements  illustrated  schematically  in  Figure 
21.  By  increasing  the  source-detector  separation  A  the  interaction  volume  is  enlarged 
while  a  minimum  value  is  enforced  for  a  wave's  paths.  According  to  our  notation,  this 


amounts  to  setting  the  lower  path-length  limit  at  .v,  =  A  and  the  upper  one  at  s2  =  A  +  S. 
Here  S  denotes  the  value  of  the  total  span  of  path-lengths  available  in  the  measurement; 
S  is  constant  in  our  experiments.  Subsequent  normalization  and  averaging  over  different 
realizations  was  performed  following  the  procedure  outlined  in  Eq.  (27). 

In  Fig.  22  we  present  the  values  of  the  normalized  scale  dependent  entropy  ha(A) 

averaged  over  ten  realizations  of  disorder  for  both  media  examined.  As  can  be  seen  for 
both  media,  when  the  interaction  volume  increases,  the  entropy  increases  as  expected 
because  in  all  realizations  a ,  8a  (s,  ga)  is  a  nonstationary  process,  and  its  fluctuations 

decrease  at  larger  5.  The  absolute  values  and  the  rate  of  increase  for  ha  ( A)  however  are 
medium  specific. 

Two  main  observations  are  in  place.  First,  we  notice  the  higher  values  of  the 
entropy  for  medium  A.  This  is  the  result  of  a  higher  number  density  Nv  of  scattering 
centers  which  determines  a  larger  number  of  possible  optical  paths  having  a  given  length 
s.  Therefore,  there  are  smaller  fluctuations  in  5 ?(s,ga)  as  discussed  before  and, 
consequently,  the  entropy  tends  toward  its  value  corresponding  to  an  infinite  number  of 
possible  trajectories  of  length  s . 


A 


-*□-*  □ 
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Figure  22  Average  normalized  entropy  /iff  (a)  for  medium  A  (blue  circles) 

and  medium  B  (red  boxes)  for  increasing  volumes  of  interaction  in  a  bistatic 
configuration  as  depicted  in  the  inset. 

The  second  observation  relates  to  the  different  rates  of  entropy  increase  as  suggested  by 
the  dashed  lines  in  Fig.  22.  This  behavior  can  be  understood  by  realizing  that  a  certain 
path-length  s  can  be  reached  through  a  different  number  m  of  scattering  events.  For 
independent  scattering,  the  joint  distribution  p(s,m )  of  such  a  process  is  Poissonian  and 


the  cumulative  probability  of  scattering  orders  up  to  M  that  contribute  to  paths  of  length 
s  is  described,  in  average,  by  a  universal  cumulative  distribution  function  p(s,M ) .  This 
cumulative  distribution  increases  fast  for  low  values  of  M  and  tends  to  saturate  for 
higher  scattering  orders.  In  one  realization  where  the  interaction  volume  is  finite,  the 
maximum  scattering  order  M  contributing  to  a  certain  5  is  essentially  determined  by  the 
number  of  available  scattering  centers.  Thus,  processes  involving  different  number 
densities  Nv  will  in  fact  experience  different  regions  of  the  cumulative  distribution 
function.  For  the  sparser  medium  B,  a  change  in  M  results  in  a  faster  increase  of  the 
corresponding  values  of  p(s,M )  and,  consequently,  a  faster  decrease  in  the  possible 
fluctuations.  Because  the  entropy  is  a  measure  of  magnitude  of  these  fluctuations,  it 
follows  that  the  medium  B  should  be  characterized  by  a  faster  rate  of  entropy  increase  as 
can  be  seen  in  Fig.  22.  As  a  result,  in  spite  of  being  described  in  average  by  the  same 
diffusion  coefficient  D ,  the  two  media  can  be  discriminated  based  on  the  value  of  Nv  . 
This  information  was  not  available  in  the  ensemble  average. 

Let  us  now  discuss  the  significance  of  these  results.  We  have  practically  demonstrated  a 
new  way  of  analyzing  the  fluctuations  of  scattered  waves  resulting  from  the  interaction 
between  coherent  fields  and  disordered  media.  In  general,  the  complexity  of  such 
deterministic  interference  can  be  reduced  only  through  ensemble  averaging  at  the 
expense  of  available  information.  We  have  demonstrated  that  analyzing  individual 
members  of  the  ensemble  of  interactions  provides  means  to  extract  information 
beyond  that  available  in  the  ensemble  average.  This  should  be  of  particular  interest  for 
remote  sensing  procudures  that  rely  on  multiple  exposures  or  extended  time-scales  in 
order  to  establish  robust  averages.  In  contrast  to  the  use  of  statistics  of  steady-state 
interactions,  the  study  of  fluctuations  in  pulse  propagation  or  interaction  with  finite 
regions  of  a  medium  open  up  additional  posibilities.  In  particular,  the  variations  from 
one  realization  of  the  interaction  to  another  can  provide  means  for  target  identifications 
or  for  fast  detection  of  structural  changes. 

The  deviation  of  an  individual  path-length  distribution  from  the  ensemble  average  is  a 
nonstationary  random  process  which  also  varies  from  one  realization  to  another.  There 
are  different  ways  to  analyze  such  a  random  function.  Here  we  have  shown  that  specific 
properties  of  the  random  medium's  morphology  can  be  evidentiated  by  using  the  scale 
dependent  entropy  associated  with  the  variance  of  path-length  fluctuations.  This  is  of 
particular  interest  in  practice  where  the  volume  of  interaction  can  be  easily  controlled 
from  macroscopic  to  microscopic  scales. 


5.  Polarization  in  randomly  scattered  fields 

Random  electromagnetic  fields  (REF)  exist  in  all  forms  and  one  common  origin  is  a 
result  of  the  interaction  of  coherent  fields  with  randomly  inhomogeneous  media.  This 
coherent  light-matter  interaction  is  a  complex  interference  process  leading  to  fields  with 
strong  fluctuations  in  intensity.  Besides  their  intensity  fluctuations,  REF  are  characterized 
by  specific  coherence  and  polarization  properties. 

A  simple  way  to  describe  the  intensity  speckles  is  to  consider  the  superposition  of  waves 
originating  from  discrete  centers  as  a  result  of  scattering.  Different  scattering  regimes 
may  vary  from  “single  scattering”  specific  to  most  surface  scattering  to  different  degrees 
of  multiple  scattering  characteristic  to  the  interaction  with  three-dimensionally  disordered 
media.  When  one  single  polarization  component  is  analyzed,  i.e.  when  the  speckle  field  is 
measured  through  a  polarizer,  the  intensity  contrast  often  reaches  unity.  This  is  the  case 
of  the  so-called  fully  developed  speckle  pattern,  a  manifestation  of  interference  between  a 
large  number  of  wavelets  with  uniformly  distributed  random  phases.  This  is  a  rather 
universal  behavior  present  in  scattering  from  a  variety  of  media  ranging  from  metallic 
rough  surfaces  to  diffusive  materials. 

However,  the  distribution  of  polarization  states  in  REF  is  much  richer  and  non-universal 
properties  are  to  be  expected.  Most  importantly,  it  is  anticipated  that  the  polarization 
properties  of  REF  corresponding  to  different  scatting  regimes  will  depend  greatly  on  the 
strength  of  the  scattering  process.  For  instance,  it  is  likely  that  when  the  wave  interaction 
is  dominated  by  single  scattering  processes,  a  fully  developed  speckle  pattern  will  occur 
but  the  REF  polarization  will  strongly  resemble  the  incident  state  of  polarization.  On  the 
other  hand,  when  the  interaction  is  subject  to  strong  multiple  scattering,  the  scattered 
field  remains  locally  fully  polarized  but  its  state  of  polarization  will  vary  from  point  to 
point. 


5.1  Scattered  fields  and  their  intensity  characteristics 

A  simple  experiment  was  designed  to  image  the  REF  at  the  surface  of  scattering  media  in 
a  backscattering  geometry.  The  samples  were  illuminated  with  a  large,  linearly  polarized 
beam  from  a  laser  with  a  wavelength  of  488nm  and  the  surface  was  imaged  onto  a  CCD. 
At  the  surface  of  the  medium,  the  speckles  are  of  the  order  of  the  wavelength  and,  in 
order  to  fully  resolve  them,  they  were  magnified  about  ninety  times  to  about  80pm  in 
size.  We  used  a  number  of  sample  media  with  varying  degrees  of  rough  surfaces  and 
volume  scattering.  Figure  23  illustrates  several  experimentally  obtained  REF 
corresponding  to  a  rough  metallic  surface  (A)  and  three  diffuse  volume  scattering  media 
characterized  by  different  transport  mean  free  paths:  a  thin  kaolin  based  diffuse  coating 
(B),  a  cellulose  membrane  (C),  and  a  polyvinylidene  fluoride  membrane  (D). 


Fig.  23.  Intensity  patterns  scattered  from  four  different  media:  a  rough  metallic  surface 
(A),  a  thin  kaolin  based  diffuse  coating  (B),  a  cellulose  membrane  (C),  and  a 
polyvinylidene  fluoride  membrane  (D).  The  insets  show  the  corresponding  log (/>(//</>)) 
vs  //</>  distributions. 


In  spite  of  their  different  origins,  the  random  fields  illustrated  in  Figure  23  are  all 
developed  speckle  patterns  as  demonstrated  by  the  probability  density  functions  of 
normalized  intensities  shown  in  the  corresponding  insets;  the  p(I/<I> )  distributions  all 
manifest  negative  exponential  decays.  In  addition,  from  the  intensity  distributions  one 
can  evaluate  the  second  order  intensity  correlations  <I(r)I(r+p)>,  i.e.  the  average  speckle 
size.  In  all  cases  shown,  the  average  speckle  size  is  approximately  the  wavelength  of 
light.  As  such,  one  can  conclude  that  both  1st  and  2nd  order  intensity  correlations  are  quite 
similar  in  all  cases  in  spite  of  the  fact  that  the  REFs  are  generated  by  scattering  on  quite 
different  random  media.  This  means  that  the  intensity  statistics  is  insufficient  to 
discriminate  between  REFs  that  are  generated  via  different  regimes  of  scattering. 

In  contrast,  it  is  known  that  the  polarization  is  more  robust  and  resistant  to  random 
fluctuations.  For  instance,  in  single  scattering  from  rough  surfaces,  the  state  of 
polarization  is  maintained.  As  the  contribution  of  multiple  scattering  increases,  the  wave 
depolarizes,  usually  over  length  scales  greater  than  the  transport  mean  free  path. 
Therefore,  it  is  expected  that  the  correlation  between  polarization  states  at  the  surface 
should  be  indicative  of  the  scattering  level  within  the  random  medium. 


5.2  Scale  dependent  degree  of  polarization 

Let  us  now  examine  the  vectorial  nature  of  the  scattered  fields  and  inspect  their 
polarization  properties.  A  rotating  quarter-wave  plate  was  inserted  in  the  imaging  path  of 
the  experimental  setup  and  a  subsequent  Fourier  analysis  provided  the  full  Stokes  vector 
in  each  pixel.  The  spatially  resolved  polarimetric  description  of  the  REF  can 
subsequently  be  analyzed  in  different  ways.  A  common  measure  of  polarization  is  the 
degree  of  polarization  (DoP) 

P(r)  =  +  S,’(r)  +  .s  (/  ./S  in ,  (30) 

S,(r)  =  E»E»-E;(r)E,(r) 

where  S2(r)  =  E*x(r)Ey(r)  -  E*y(r)Ex(r)  .  (31) 

S3(r)  =  i(E'x(r)Ey(r)-E;(r)Ex(rj) 

As  can  be  seen,  the  DoP  is  actually  a  4th-order  field  correlator  describing  the  properties  of 
the  electric  field  at  point  r.  Practically,  one  has  access  to  an  ensemble  of  different 
polarization  states  collected  over  a  certain  area  and  a  scale  dependent  effective  DoP  can 
be  estimated  as 


P  (r)  =  If  Stdr  +  j  S'dr  +  j  S'dr  /\Idr  ■  @2) 

V  A  A  A  A 

The  scale  dependent  DoP  in  Eq.  (32)  approaches  zero  when  the  ensemble  of  polarization 

states  are  randomly  distributed  over  the  Poincare  sphere  as  illustrated  in 

Fig.  14.  Of  course,  a  strong  multiple  scattering  interaction  leads  to  an  overall  depolarized 
scattered  field  when  the  DoP  is  evaluated  over  a  large  scale.  We  emphasize  that  this 
depolarization  only  occurs  in  the  global  sense,  as  the  resulting  scattered  light  forms 
individually  fully  polarized  speckles.  The  depolarization  measurement  can  be  viewed  as  a 
“center  of  mass”  estimation,  where  the  polarization  of  the  average  state  lies  within  the 
sphere  as  seen  in  Fig  4b. 


Fig.  14.  a)  Distribution  of  polarization  states  on  the  Poincare  sphere,  b)  Both  the 
average  state  of  polarization  and  the  degree  of  polarization  are  ensemble  properties  of 
the  distribution  of  polarization  states. 


Results  of  calculations  of  scale-dependent  DoP  are  presented  in  Fig.  25.  As  we  increase 
the  size  of  the  area  over  which  the  estimation  is  performed,  the  DoP  value  eventually 
reaches  saturation.  The  error  bars  in  Fig  25  were  calculated  using  one  hundred  starting 
points  for  each  size  of  analysis  area.  Their  values  reflect  the  fact  that  even  if  each  speckle 

is  fully  polarized,  PA  (r)  varies  significantly  when  only  examining  a  small  number  of 
speckles. 

As  evident  from  Fig.  25,  there  are  two  factors  associated  with  the  DoP  scale  dependence: 
the  saturation  level  and  the  associated  decay  length.  It  is  quite  clear  that  the  saturation 
levels  can  be  quite  different.  For  the  samples  illustrated  here,  the  DoP  saturation  values 
are  0.97,  0.49,  0.32,  and  0.32  +/-  0.02  for  samples  A,  B,  C,  and  D,  respectively.  Different 
saturation  levels  indicate  different  levels  of  global  depolarization  due  to  scattering:  the 
higher  the  level  of  multiple  scattering  the  lower  the  value  of  the  corresponding  degree  of 
global  polarization  of  the  scattered  field.  However,  the  two  membrane  samples  (C  and  D) 
show  very  similar  DoP  saturation  values  yet  they  have  completely  different  material 
structures.  In  this  case  the  DoP  value  is  insufficient  to  discriminate  between  different 
levels  of  multiple  scattering. 


Fig.  25.  The  DoP  calculated  over  an  increasing  integration  area  averaged  over  100  integration 
centers.  The  error  bars  indicate  the  DoP  standard  deviation  calculated  for  each  integration 
area.  Note  that  the  y-scale  in  A  is  0.95-1.0  while  it  is  0-1.0  for  panels  B,  C,  and  D. 

The  DoP  decay  length  represents  the  spatial  scale  at  which  saturation  is  reached. 
Although  the  media  examined  here  have  different  saturations,  they  all  have  similar  decay 


lengths  measuring  about  15  intensity  speckles.  It  appears  that  this  decay  length  is  an 
associated  ensemble  quantity  that  ignores  any  underlying  material  discrepancies. 

To  conclude,  the  DoP  is  a  basic  yet  incomplete  way  to  describe  the  polarization 
properties  of  REF.  As  mentioned  before,  the  DoP  is  evaluated  based  on  4th-order  field 
correlations  representative  for  the  properties  of  a  field  at  one  spatial  location.  It  is  also  an 
ensemble  quantity,  averaged  over  many  spatial  locations,  and  therefore  it  loses  any 
information  about  spatial  distribution  of  polarization  states,  i.e.  the  shape  of  the 
distribution  of  states  on  the  Poincare  sphere.  For  instance,  there  are  many  different  REF 
that  can  be  characterized  by  the  same  value  of  the  global  DoP.  In  the  following  we  will 
examine  a  more  refined,  two-point  descriptor  of  the  polarization  properties  of  a  random 
field. 


5.3  Two  point  polarization  correlations 

Another  possibility  to  evaluate  higher-order  field  correlations  in  REF  is  to  quantify  the 
variability  between  the  members  of  the  ensemble  of  measurement  points  and  a  chosen 
reference.  The  measure  for  this  is  the  complex  degree  of  mutual  polarization  (CDMP)4. 
The  magnitude  of  the  CDMP  measures  the  polarization  similarity  between  two  points,  r, 
and  a  reference  r0,  and  under  the  assumption  of  a  fully  correlated  and  locally  fully 
polarized  field,  it  is  defined  as 


V2(^,r0)  = 


[K^em  +  eI^e^))2 


Ex(0\2 +\Ey(r0)\  \Ex(rtf +\Ey(r0) 


(33) 


As  such,  CDMP  measures  the  distance  between  two  points  (states  of  polarization)  on  the 
Poincare  sphere;  orthogonal  states,  opposite  to  each  other  on  the  Poincare  sphere,  have  an 
associated  CDMP  value  of  zero.  The  CDMP  reflects  the  shape  of  the  distribution  of 
states  on  the  Poincare  sphere  in  contrast  to  the  DoP,  which  is  a  measure  of  the  location  of 
the  center  of  mass  of  the  distribution.  CDMP  quantifies  the  spatial  distribution  of  states 
by  comparing  the  Stokes  vectors  in  each  point  to  a  common  reference.  Recently,  similar 
estimations  of  REF  properties  permitted  to  detect  local  non-stationarities  such  as  the 
presence  of  a  weak  localization  phenomenon. 

It  is  important  to  note  that  CDMP  is  not  an  ensemble  quantity  and  it  can  be 
calculated  while  maintaining  spatial  information.  In  Fig  26  we  show  the  spatial 
distribution  of  the  REFs  in  Fig.  23  but  this  time  encoded  in  the  CDMP  values  calculated 
with  respect  to  the  constant  state  of  polarization  of  the  incident  field. 


4  J.  Ellis  and  A.  Dogariu,  "Complex  degree  of  mutual  polarization,"  Opt.  Lett.  29,  536-538  (2004). 
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Fig.  26.  CDMP  maps  of  the  REF  in  Fig.  1  calculated  with  respect  to  incident  state  of 
polarization.  Insets  show  binary  images  of  corresponding  maps  thresholded  at 


CDMP=0.5. 


As  can  be  seen,  the  CDMP  maps  provide  information  about  the  similarity  of  polarization 
states  in  the  REF.  The  rough  metallic  surface,  sample  (A),  exhibits  strong  single 
scattering  and  as  a  result  most  of  the  scattered  field  is  in  the  same  state  of  polarization  as 
the  incident  one;  the  CDMP  is  almost  uniformly  unity  as  a  result  of  the  strong  spatial 
correlation  of  the  polarization.  As  the  level  of  multiple  scattering  increases,  the  spatial 
correlation  of  the  CDMP  maps  reduces  as  seen  for  samples  B,  C,  and  D. 

To  assess  quantitatively  the  differences  between  the  CDMP  maps,  one  can  examine  the 
probability  density  functions  (PDF)  corresponding  to  the  distributions  of  CDMP  values 
across  the  image.  As  can  be  seen  in  Fig.  27,  different  distributions  are  found  for  the 
CDMP  maps  shown  in  Fig  26.  Clearly,  the  PDFs  evolve  from  being  peaked  around  unity 
in  the  case  of  single  scattering  (rough  surface  A)  to  a  more  uniform  distribution 
corresponding  to  the  higher  order  multiple  scattering  processes  in  volume  scattering 
media.  An  interesting  observation  can  be  made  in  regard  to  samples  C  and  D.  Even 
though  the  REF  scattered  by  these  media  have  almost  the  same  global  degree  of 
polarization,  the  PDFs  of  their  corresponding  CDMP  maps  are  noticeably  different.  This 
can  be  interpreted  as  different  coverage  on  the  Poincare  sphere  corresponding  to  similar 
centers  of  mass.  Being  a  local  4th-order  field  correlation,  CDMP  preserves  the  spatial 
information  and  distinguishes  between  different  underlying  field  distributions  that  have 
similar  average  properties. 


Fig.  27  Probability  density  functions  for  the  CDMP  distribution  maps  shown  in  Fig.  26. 

We  note  that  for  a  fully  diffusive  process  of  wave  interaction  all  polarization  states  are 
equally  probable  and  a  uniform  PDF  is  to  be  expected.  The  slight  increase  of  the 
probability  density  towards  higher  CDMP  values  means  that  the  diffusive  scattering 
behavior  has  not  been  reached.  In  the  present  case  this  is  a  consequence  of  the 
backscattering  geometry  where  the  low-order  scattering  events  are  always  dominant. 

Let  us  now  turn  our  attention  to  the  spatial  characterization  of  the  CDMP  maps  in  Fig.  26. 
As  can  be  seen,  these  maps  reveal  areas  of  more  uniform  polarization.  One  can  interpret 
theses  areas  of  spatial  variation  of  the  CDMP  values  as  determining  some  sort  of 
polarization  or  CDMP  “speckles”.  This  spatial  correlation  over  the  CDMP  maps  is  more 
obvious  in  the  binary  images  included  as  insets  in  Fig.  26.  The  binary  distributions  were 
obtained  by  thresholding  all  the  points  having  a  CDMP  value  greater  than  0.5.  The  field 
in  A  is  clearly  dominated  by  single  scattering,  where  the  incident  state  of  polarization  is 
largely  maintained.  In  the  binary  images  corresponding  to  the  other  samples,  the  different 
extent  of  clearly  defined  areas  with  similar  polarization  is  evident,  which  can  be 
understood  as  different  levels  of  spatial  correlations  between  the  CDMP  values.  This  is 
akin  to  defining  an  average  size  of  the  CDMP  “speckle”. 

When  considering  the  distributions  in  the  binary  map  shown  in  Fig.  26B,  we  found  that 
the  average  size  of  correlated  areas  is  about  5.3  times  larger  than  the  average  size  of  an 
intensity  speckle.  The  size  of  the  correlation  areas  clearly  decreases  for  samples  C  and  D 
where  it  is  3.3  and  3.0  intensity  speckles,  respectively.  Of  course,  the  number  of  clusters 
and  their  size  depend  upon  the  threshold  level.  As  the  threshold  decreases,  the  total 
number  of  points  that  lay  above  this  level  increases  and  the  size  of  these  correlated  areas 
increases  as  shown  in  Fig.  28.  However,  the  areas  increase  at  different  rates  as  clearly 
seen  in  Fig.  28,  indicating  different  structures  in  the  CDMP  maps  and,  moreover,  a 
distinctive  dependence  on  the  level  of  scattering  in  each  sample.  We  should  reiterate  that 


both  the  distribution  of  CDMP  values  in  Fig.  27  and  the  extent  of  the  two-dimensional 
CDMP  correlations  in  Fig.  6  are  obtained  from  one  single  realization  of  the  random  fields 
shown  in  Fig.  23.  Thus,  this  analysis  constitutes  the  point  and  the  point-pair 
characterization  of  the  polarization  properties  of  the  scattered  electromagnetic  field  in  a 
manner  similar  to  the  intensity  distribution  p(I)  and  the  average  intensity-intensity 
correlation  <I(r)I(r+R)>  in  a  scalar  speckle  pattern. 
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Fig.  28  Size  of  CDMP  speckle  in  units  of  intensity  speckles  as  a  function  of  CDMP 
threshold  level.  Different  decay  rates  are  evident  for  the  four  samples  examined. 


Finally,  we  would  like  to  briefly  comment  on  the  information  contained  in  these  REF 
properties.  It  is  known  that  a  number  of  polarization  memory  effects  are  present  at 
different  levels  of  scattering.  There  is  an  intimate  dependence  between  the  medium’s 
structure  and  the  polarimetric  properties  of  the  scattered  field  and,  therefore,  one  can 
anticipate  that  the  distribution  of  polarization  states  and  their  spatial  correlation  in  a  REF 
should  reflect  some  of  the  morphological  properties  of  the  scattering  media.  Let  us 
consider  again  the  two  samples  that,  in  average,  depolarize  the  light  at  essentially  the 
same  level,  DoP  ~  0.32 ,  yet  their  structural  morphology  is  quite  different.  From  Figs. 
27  and  28,  one  can  clearly  see  that  both  the  PDF  of  the  corresponding  CDMP  maps  and 
the  sizes  of  the  CDMP  speckles  are  different  for  these  samples.  This  is  because  the 
structural  differences  lead  to  different  scattering  strengths  in  these  two  media.  To  asses 
these  differences  we  performed  typical  ensemble  average  measurements  of  enhanced 
back  scattering  (EBS).  These  measurements  yielded  different  values  of  the  transport 
mean  free  path:  8  and  7  pm  for  media  C  and  D,  respectively,  as  estimated  from  the  full 
width  at  half  maximum  of  the  enhancement  peaks. 


Being  a  measure  of  polarization  similarity  at  different  spatial  locations,  the  size  of  the 
CDMP  speckle  reflects  the  extent  of  the  interaction  volume  necessary  for  the  wave  to 
depolarize,  or,  in  other  words,  to  lose  memory  of  its  initial  polarization  state.  In  a  specific 


geometry,  the  magnitude  of  this  characteristics  length  scale  depends  on  the  number  of 
transport  mean  free  paths.  Therefore,  scattering  media  characterized  by  small  values  of  1* 
are  also  expected  to  generate,  at  their  surface,  scattered  fields  with  smaller  value  of  the 
CDMP  speckles.  This  is  exactly  what  our  experiments  show;  the  lowest  value  of  the 
CDMP  speckle  corresponds  to  the  strongest  scattering  in  sample  D.  Remarkably,  one 
single  realization  of  the  scattering  process  is  sufficient  to  provide  information  similar  to 
that  acquired  through  an  ensemble  average  measurement. 

In  conclusion,  we  have  found  that  the  spatial  variability  of  the  vectorial  properties  can  be 
markedly  different  even  when  the  random  fields  have  similar  global  properties.  The  point 
and  point-pair  correlations  of  the  complex  degree  of  mutual  polarization  provide  means 
to  identify  the  origins  of  scattered  fields.  We  demonstrated  that  the  extent  of  these  spatial 
correlations  is  determined  by  the  magnitude  of  the  transport  mean  free  path 
characterizing  the  scattering  process.  Spatially  resolved  measurements  of  the  polarization 
properties  in  one  realization  of  the  scattered  field  allow  recovering  information  otherwise 
available  only  through  ensemble  averages. 


6.  Unique  polarization  effects  in  scattering 

One  of  the  important  features  of  active  systems  is  that  one  can  selectively  probe  different 
object  characteristics  by  controlling  the  incident  illumination.  A  possibility  of  interest  is 
to  structure  the  illuminating  beams  in  terms  of  its  polarization,  intensity,  coherence,  or 
phase.  It  is  expected  that  the  structure  of  the  illuminating  field  can  greatly  influence  the 
information  probed.  For  example,  by  illuminating  a  diffusive  media  with  phase  vortices 
of  differing  charge,  different  path-lengths  are  selectively  probed  in  the  enhancement 
peak5. 

In  order  to  take  advantage  of  such  unique  features  one  has  to  understand  the  details  of 
propagation  and  scattering  of  polarized  radiation.  In  particular,  there  is  a  growing  interest 
in  using  the  specific  properties  of  circularly  polarized  light  (robustness  in  propagation, 
possibility  of  digital  processing  and  others).  Here  we  report  on  novel  characteristics  of 
scattering  of  circularly  polarized  light.  These  novel  effects  can  be  speculated  for  the 
design  of  new  types  of  radiation  sensors  and  for  interpreting  specific  scattering 
signatures.  In  addition,  this  research  may  also  lead  to  new  possibilities  for  optical 
manipulation  of  matter. 

In  elastic  scattering  from  a  non  absorbing  spherical  particle,  two  parameters  of  the 
electromagnetic  field  are  conserved:  the  energy  and  the  angular  momentum  component 
along  the  propagation  direction.  Linear  momentum  as  a  whole  is  also  conserved,  of 
course,  but  some  of  it  is  transferred  to  the  particle  leading  to  radiation  pressure.  The 
conservation  of  energy  (which  is  a  scalar  quantity)  leads  to  a  normalization  condition  for 
the  integrated  energy  flux  density,  which  is  further  used  in  defining  the  scattering  cross 
section.  Similarly,  the  conservation  of  angular  momentum  should  be  related  to  the 
angular  momentum  flux  density.  The  continuity  conditions  for  the  angular  momentum 
density  can  be  described  by  three  equations  (one  for  each  component  of  the  vector)  or 
one  equation  for  a  tensor. 

We  are  interested  in  calculating  the  angular  moment  flux  density  of  the 
electromagnetic  field  which  results  from  scattering  of  a  circularly  polarized  wave  from  a 
non-absorbing  spherical  particle.  It  was  demonstrated  that  scattering  of  circularly 
polarized  wave  does  not  exert  torque  on  the  particle  and  that  transfer  of  angular 
momentum  from  the  field  to  the  particles  is  mediated  only  by  absorption.  Therefore,  the 
angular  momentum  of  the  field  is  preserved  and  it  should  be  interesting  to  know  how  the 
flux  density  of  the  angular  momentum  is  distributed  between  the  spin  or  polarization  term 
which  will  be  henceforth  designated  by  s  and  the  orbital  angular  momentum  (OAM), 
which  will  be  denoted  by  1 .  The  total  angular  momentum  flux  density,  j,  is  given  by  the 
sum  of  these  two  terms. 


5  Field  correlations  under  vortex  illumination  have  been  found  to  act  as  a  type  of  path-length  filter. 
Furthermore,  modification  of  the  illuminating  structure  allows  for  a  deterministic  change  in  the  measured 
field  correlations,  allowing  for  a  novel  and  simplified  means  of  measuring  structural  properties  of  the 
medium,  cf.  C.  Schwartz  and  A.  Dogariu,  “Coherent  backscattering  of  singular  beams,”  Opt.  Commun  258 
(2006) 


6.1  Scattering  effects  on  spin  angular  momentum 

Let  us  consider  an  incident  plane  wave  which  is  monochromatic  (angular  frequency  co), 
left  circularly  polarized,  has  an  amplitude  E0  and  propagates  in  the  direction  z  (the  'hat' 

denotes  a  unit  vector).  This  is  the  simplest  example  of  a  paraxial  wave  which  has  the 
general  form 


E(x,  y,  z )  =  exp  (ikz)  F(x,  y,  z)  (34) 

(where  F(x,  y,z)  is  a  slowly  varying  spatial  envelope).  For  such  a  wave  one  can  employ 
expressions  for  angular  momentum  flux  density  of  a  paraxial  wave  and  write  the  z 
component  of  the  angular  momentum  flux  density  can  be  written  as 
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The  first  term  relates  to  the  transverse  distribution  of  the  field  and  is  the  OAM  term.  The 
second  term  reflects  the  angular  momentum  carried  by  circular  polarization  is  zero  for 
linear  polarization  and.  For  a  class  of  paraxial  beams  for  which  the  transverse  field 
distribution  can  be  written  as 


and  recognizing  that 
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jz(r’0)=^T-(m  +  (7) \u(rf  ={m  +  o)—.  (37) 

2  co  co 

In  Eq.  (37)  cr  is  the  wave  helicity:  +1  for  left  circular  polarization,  -1  for  right  circular 
polarization  and  0  for  any  linear  polarization.  These  beams  are  called  vortex  beams  and 
they  exhibit  an  axial  singularity  since  the  phase  is  indeterminate  on  the  beam  axis.  Due  to 
destructive  interference,  this  phase  singularity  leads  to  an  intensity  null.  The  integer 
factor  m  is  called  the  topological  charge  of  the  singularity  and  it  is  an  invariant  of  the 
propagation.  This  follows  from  the  fact  that  m  relates  directly  to  the  orbital  angular 
momentum  and,  under  the  paraxial  approximation,  each  term  of  the  angular  momentum  is 
conserved  independently. 

The  flux  density  of  the  angular  momentum  of  the  incident  plane  wave  which  is  carried  by 
the  spin  term  is  equal  to  I0/co  where  /0  =  s0cEq/2  ( e0  being  the  permittivity  of  vacuum 

and  c  the  velocity  of  light).  This  corresponds  to  a  quantum  description  in  which  each 
photon  in  the  incident  field  carries  a  spin  angular  momentum  of  +  y .  The  rate  at  which 


the  angular  momentum  is  removed  from  the  incident  field  can  be  derived  from  the  the 
scattered  power  and  is  given  by  <jscIJ col  .  This  is  the  source  term  for  the  angular 

momentum  of  the  scattered  field  and  it  should  be  recovered  by  integrating  over  all 
directions  the  angular  momentum  flux  density  of  the  scattered  field. 

From  general  scattering  theory  it  is  known  that  the  scattered  electric  field  in  the  far  zone 
can  be  written  as 
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where  the  matrix  S  is  the  scattering  matrix  in  the  circular  polarization  basis,  the  L  and 
R  subscripts  designate  the  left  and  right  circular  polarizations.  In  Eq  (38),  k  is  the  wave 
number  and  r  is  the  distance  from  the  center  of  the  scattering  particle  which  is  both  the 
origin  of  the  coordinates  frame  and  the  reference  point  for  the  angular  momentum 
calculations.  The  scattering  matrix  in  the  circular  base  can  be  related  to  the  more  usual 
amplitude  scattering  matrix  given  in  terms  of  the  parallel  and  perpendicular  electric  field 
components  with  respect  to  the  scattering  plane.  For  a  spherical  particle,  this  scattering 
matrix  is  diagonal  with  elements  S2{0),  5, (d)  and  the  relation  is 
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The  scattered  field  can  now  be  written  as 


E (6, <f>)  =  exP^r)  eo  exp  (/^J 


(40) 


and  expanding  further  the  left  and  right  circular  unit  vectors  in  terms  of  the  locally 
transverse  unit  vectors 
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one  finally  obtains 
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where  the  notation 

-se (6>)  =  (<?)+  (0)VV2  and  S#(e)=i[5It(f?)-5„(^,/2  has  been  used. 


Having  found  E,  we  can  now  proceed  to  calculate  the  angular  momentum  flux  density 
through  a  radially  oriented  infinitesimal  area  in  the  radiation  zone.  When  the  wave  is 
approximated  locally  as  a  plane  wave  the  spin  term  can  be  found  to  be 
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This  result  is  physically  reasonable:  the  spin  is  simply  the  difference  of  intensities  of  the 
two  radially  outgoing  orthogonal  circular  polarization  components  divided  by  the  angular 
frequency  of  radiation.  The  common  phase  term  exp(z^)  does  not  play  a  roll  in  this 
calculation. 


Due  to  axial  symmetry,  all  the  components  of  s  average  to  zero  when  an  integration  over 
the  angles  is  performed  except  for  the  z  component  which  is 
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For  a  Rayleigh  scatterer,  one  can  immediately  find  that  the  expression  in  Eq.  (44)  reduces 
to 


(o)  =  \  °sc  cos2  (<?)•  (45) 
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In  order  to  calculate  the  total  scattered  flux,  the  spin  flux  density  must  be  integrated  over 
a  sphere  of  radius  r  to  obtain 
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The  bar  superscript  in  Eq.  (46)  indicates  integration  over  a  sphere  of  arbitrary  radius,  in 
the  far  field.  Notably,  one  can  see  that  only  half  of  the  angular  momentum  flux  removed 


from  the  incident  wave  is  contained  in  the  spin  term.  Of  course,  this  conclusion  is  similar 
to  the  results  obtained  for  the  case  of  a  field  radiated  by  a  rotating  dipole.  A  similar 
treatment  can  be  developed  for  a  the  more  complex  case  of  a  Mie  scatterer. 


6.2  Scattering  effects  on  orbital  angular  momentum 

Let  us  now  turn  our  attention  to  the  more  complex  issue  of  the  OAM  term  which  is 
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In  order  to  evaluate  this  expression  we  employ  the  spherical  coordinates  form  of  the 
gradient  operator  and  recall  that  the  unit  vectors  have  to  be  differentiated  as  well;  for 
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phase  term,  exp(i^),  is  important  in  this  case  and  should  be  considered  in  the 
differentiation.  Accounting  only  for  the  terms  that  contribute  to  the  z  component,  we 
obtain  after  some  algebra  that 


This  result,  together  with  the  one  expressed  in  Eq.  (46)  indicate  that  s,(0)  and  L (d)  sum 
up  to  an  expression  which  is  proportional  to  the  scattered  intensity.  For  each  6,  the  total 
angular  momentum  flux  density  is  simply  the  intensity  divided  by  the  angular  frequency 
of  the  radiation.  Integrating  over  a  spherical  surface  in  the  far  field  leads  indeed  to  a 
manifestation  of  the  conservation  of  angular  momentum  flux: 

jz=sz+lz=C7sc—[z-  (49) 

CO 


Moreover,  we  note  again  that  at  each  6 ,  the  ratio  j_  (0)/l(@)  is  constant.  In  the  quantum 
description  one  may  say  that  a  photon  scattered  in  any  direction  carries  the  same  angular 
momentum  as  an  incident  photon,  but,  in  different  directions  it  is  distributed  differently 
between  the  spin  and  OAM  terms.  It  is  now  evident  that  in  the  forward  and  backward 
directions  the  angular  momentum  is  carried  only  by  the  polarization.  In  other  words  the 
helicity  is  fully  preserved  in  the  forward  direction  and  it  is  fully  reversed  in  the  backward 
direction.  For  some  angles,  which  for  a  Rayleigh  scatterer  is  only  90 0 ,  on  the  other  hand, 
the  scattered  light  is  linearly  polarized  and  therefore  does  not  carry  any  spin  angular 
momentum.  At  these  angles  the  angular  momentum  is  carried  only  by  OAM  term. 
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Figure  29.  Three-dimensional  angular  distributions  of  the  normalized  spin  term 
(angular  momentum  content  of  the  spin  term  divided  by  the  intensity)  for 
several  size  parameters.  The  relative  index  of  refraction  was  1.09  .  Left 
circularly  polarized  light  is  incident  along  the  z  axis  (the  direction  is  indicated 
by  the  red  arrow).  The  color  bars  represents  the  spin  and  a  complementary 
illustration  can  be  obtained  for  the  OAM 

Figure  29  illustrates  the  three  dimensional  distribution  of  the  spin  term  (relative  to  the 
scattered  intensity)  plotted  as  a  function  of  scattering  direction  for  several  values  of  the 
size  parameter  x  =  27ml A, ,  where  a  is  the  radius  of  the  particle  and  A  is  the  wavelength. 
The  orbital  angular  momentum  distribution  is  the  complementary  one  (as  it  will  be  given 
by  one  minus  the  spin).  It  is  interesting  to  note  that  for  small  size  parameters  the  changes 
in  the  normalized  distribution  are  small.  For  large  particles  on  the  other  hand,  in  which 
Mie  resonances  are  dominant,  one  can  observe  that  there  are  directions  in  which  the  spin 
term  is  negative-  meaning  that  the  normalized  orbital  term  is  larger  then  1.  The 
transition  to  such  behavior  happens  for  particles  with  a  size  parameter  of  about  n  (i.e.  a 
particle  diameter  which  is  approximately  equal  to  the  wavelength). 

In  order  to  illustrate  the  distinction  between  spin  and  orbital  angular  momenta  one  may 
use  the  following  "gedanken"  experiment.  Let  us  consider  a  small  dielectric  sphere  which 
is  slightly  absorbing,  and  is  placed  in  the  far  field  of  a  Rayleigh  scatterer.  In  the  exact 
forward  direction  the  particle  will  rotate  about  itself  due  to  the  absorption  of  circularly 
polarized  light.  At  of  90 0  with  respect  to  the  direction  of  incidence,  on  the  other  hand,  the 
test  particle  will  "orbit"  the  scattering  particle  due  to  the  phase  gradient  in  the  scattered 
field.  We  would  like  to  emphasize  that  this  effect  is  not  due  to  absorption  and  that,  of 
course,  the  test  particle  will  experience  a  radial  force  due  to  radiation  pressure. 


Figure  30.  Qualitative  illustration  of  the  Poynting  vector  direction  for 
radiation  scattered  by  a  Rayleigh  scatterer  illuminated  by  a  circularly 
polarized  light  which  is  perpendicular  to  the  plane  of  the  figure  as  indicated. 

The  spiraling  of  the  Poynting  vector  results  in  an  apparent  angular  shift  of  the 
light  at  the  far  field. 

Finally,  we  note  that  the  orbital  angular  momentum  of  scattered  light  can  also  be 
interpreted  as  the  result  of  a  slight  direction  dependant  shift  of  the  apparent  origin  of  the 
scattered  waves,  a  shift  which  introduces  an  "impact  parameter"  of  the  order  of  XII. n. 
The  origin  of  this  shift  can  be  identified  in  the  classical  electromagnetic  theory  where  it  is 
known  that  the  Poynting  vector  of  radiation  from  a  rotating  dipole  (equivalent  to  a 
Rayleigh  scatterer  illuminated  by  a  circularly  polarized  light)  spirals  in  the  near  field.  As 
a  result  the  far  field  radiation  seems  to  be  emerging  not  from  the  center  of  the  dipole  but 
from  a  shifted  position  as  illustrated  qualitatively  in  Figure  30.  The  angle  a  is 
approximated  by  AJ2ur  and  the  angular  momentum  at  the  point  r  can  then  be  calculated 
to  be  kra  =  1 .  This  means  of  course  that  the  orbital  angular  momentum  carries  all  the 
angular  momentum.  In  our  gedanken  experiment,  the  slight  tilt  of  the  k-vector  induces  a 
rotating  motion  of  the  test  particle  placed  at  r .  This  is  the  mechanism  which  couples 
angular  momentum  from  the  electromagnetic  field  to  the  medium. 


6.3  Spin  transport  through  dielectric  interfaces  -  spin  Hall  effects  of 
light 

We  conducted  an  experiment  to  demonstrate  one  important  consequence  of  this  spin-orbit 
interaction.  As  we  have  seen,  in  addition  to  energy  and  linear  momentum,  a  light  wave 
carries  angular  momentum  distributed  between  orbital  E'trxVlEl  and  spin 
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macroscopic  energy  flow  with  respect  to  an  arbitrary  reference  point  and  does  not  depend 
on  the  wave’s  polarization.  The  spin  angular  momentum  (SAM)  on  the  other  hand, 
relates  to  the  phase  between  orthogonal  field  components  and  is  completely  determined 
by  the  state  of  polarization.  Because  the  total  angular  momentum  is  conserved  along 
propagation,  when  a  wave  encounters  a  refractive  index  inhomogeneity  (including 
interfaces  between  otherwise  homogeneous  media),  a  momentum  transfer  appears 
between  the  orbital  and  spin  terms.  This  angular  momentum  exchange,  or,  in  other 
words,  the  spin-orbit  interaction  (SOI),  explains  a  number  of  polarization  effects.  But 
perhaps  most  attractive  is  that  this  description  of  electromagnetic  wave  propagation 


provides  an  analogy  with  electron  transport  and  related  transport  phenomena  such  as  the 
spin  Hall  effect  in  semiconductors. 
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Figure  31.  (a)  A  beam  of  light  incident  on  a  planar  refractive  index  gradient,  (b) 
beam  shifts  observed  in  far-field  corresponding  to  incidence  in  pure  states  of 
polarization  (right  circular,  left  circular),  and  (c)  the  observable  result 
corresponding  to  an  incident  beam  in  a  mixed  state,  (d)  A  plane  wave  incident  on 
a  planar  refractive  index  gradient,  (e)  the  far-field  shifts  corresponding  to 
incidence  in  pure  states  of  polarization,  and  (f)  the  result  of  an  incident  plane  wave 
in  a  mixed  state  of  polarization,  (g)  A  plane  wave  scattering  from  a  sphere,  (h)  the 
transversal  shifts  in  the  perceived  sphere  center  as  observed  in  far-field  for  pure 
incident  states,  and  (i)  the  result  of  scattering  for  incident  wave  in  a  mixed  state. 


There  is  a  host  of  Hall-like  effects  in  electronic  systems,  all  dealing  with  transport 
phenomena  resulting  in  some  displacement  current.  Recently,  it  was  suggested  and 
demonstrated  that  an  equivalent  spin  Hall  effect  of  light  (SHEL)  exists  and  can  be 
measured  for  a  beam  impinging  on  a  dielectric  interface  (Fig.  31a).  In  SHEL,  the  electron 
spin  is  replaced  by  the  wave’s  polarization,  and  the  role  of  an  applied  electric  field 
generating  the  electronic  current  is  taken  by  the  refractive  index  gradient.  The  effect  is 
observed  as  a  displacement  of  the  beams  carrying  spin,  i.e.  a  transversal  shift  in  the 
perceived  location  of  the  interaction  volume.  For  different  incident  spins,  the  shift  is  in 


opposite  directions,  which  is  analogous  to  the  induced  electron  spin  flux  perpendicular  to 
the  initial  electronic  current.  When  the  incident  beam  is  in  a  pure  state  (circularly 
polarized),  a  shift  in  the  beam’s  center  of  mass  may  be  observed  as  illustrated  in  Fig.  31b. 
Exciting  with  a  mixed  state,  results  in  a  separation  of  spin,  where  the  region  of  overlap  is 
still  in  a  mixed  state  of  polarization  as  suggested  in  Fig  31c. 

We  will  now  examine  the  case  where  a  homogeneous  plane  wave  is  incident  on  the 
planar  dielectric  interface  as  depicted  in  Fig.  3 Id.  Due  to  the  infinite  extent  of  the 
interaction,  the  transverse  shift  cannot  be  observed  (Fig.  31e),  as  the  refracted  field  is  still 
an  infinite  plane  wave  and  no  reference  point  can  be  identified.  When  the  incident  field  is 
not  in  a  pure  state,  for  example,  linearly  polarized,  in  spite  of  the  transversal  spin  fluxes, 
there  is  no  net  photon  current  and  therefore  there  will  be  no  observable  effect  (Fig.  3 If). 

Let  us  now  inspect  cases  of  higher  dimensionality.  If  the  refractive  index  gradient  is 
spherical  and  the  excitation  is  a  plane  wave  as  illustrated  in  Fig.  8g,  the  conservation  of 
angular  momentum  results  in  transversal  spin  flows  tangent  to  the  spherical  surface.  In 
this  case,  performing  sequential  excitations  with  pure  states  and  using  a  detection  scheme 
sensitive  to  the  local  direction  of  energy  flow,  the  transversal  shifts  can  be  observed  as 
we  will  show  later.  The  presence  of  this  shift  breaks  the  symmetry  of  the  field 
distribution  relative  to  any  plane  that  contains  the  propagation  vector  and  the  center  of  the 
sphere.  In  other  words,  the  result  of  the  interaction  between  circularly  polarized  light  and 
a  sphere  depends  on  the  incident  spin:  the  sphere  distinguishes  between  left  and  right. 
For  excitation  with  a  mixed  state,  the  transversal  shift  that  occurs  cannot  be  directly 
observed  because,  again,  there  is  no  net  transverse  photon  flux  as  in  the  case  of  plane 
wave  impinging  on  a  plane  interface.  It  should  be  noted  that  the  similar  considerations 
may  be  applied  in  a  cylindrical  geometry.  These  circumstances  are  similar  to  the  Corbino 
disk  geometry  in  electronic  systems,  where  counter-propagating  spin  currents  are 
generated  but  no  net  electron  fluxes  can  be  detected. 

This  spherical  geometry  warrants  a  closer  examination.  It  is  now  rather  well  understood 
that  when  circularly  polarized  light  scatters  from  a  dipole,  the  far-field  Poynting  vector  S 
makes  an  angle  with  the  line  of  direct  sight  (radial  vector  p).  This  determines  a  shift  in 
the  perceived  position  of  the  dipole.  The  shift  is  perpendicular  to  the  plane  of  scattering, 
and  has  opposite  signs  for  excitations  with  different  spins  as  seen  in  Fig  32.  Its  value 

^oie  =  4  /  2^)[sin  0 /( 1  -  l/2sin2  0\  (50) 

depends  on  the  polar  angles  of  scattering  6  and  reaches  a  maximum  in  the  plane  where  a 
full  transformation  from  SAM  to  OAM  occurs:  the  plane  perpendicular  to  the  direction  of 
propagation  of  the  incident  wave  {6  =  ;r/2).  Note  that  this  shift  does  not  depend  on  the 
optical  properties  of  the  scatterer.  For  a  sphere  of  finite  size,  in  the  far  field,  the  shift  in 
perceived  location  is  A  =  limp(|sJ/|s  |)>  where  S_  and  S  are  the  components  of  S 

perpendicular  and  parallel  to  p. 

In  general  the  transformation  of  SAM  to  OAM  is  not  complete.  However,  for  certain  sets 
of  parameters  (sphere  radius,  refractive  index)  there  are  one  or  several  angles  6  where 
transformation  is  total  and  the  scattered  light  is  linearly  polarized.  The  value  of  A 
corresponding  to  these  angles  attains  local  maxima  and  can  reach  tens  of  wavelengths  in 


magnitude.  This  modification  in  the  perceived  location  of  the  sphere  as  shown  in  Fig.  32 
is  the  manifestation  of  the  spin  Hall  effect  in  scattering  from  a  finite  object. 


Figure  32.  Poynting  vector  field  lines  for  radiation  scattered  by  a  Rayleigh 
scatterer  illuminated  by  left  and  right  circularly  polarized  light.  The  spiraling 
of  the  Poynting  vector  results  in  the  far-field  angular  shift  suggested  in 
Figure  30. 

There  is  one  feature  in  the  scattering  of  circularly  polarized  light  from  a  sphere  that  has 
no  analog  in  previous  descriptions  of  spin  transport  effects:  Poynting  vector  circulation 
around  the  sphere  in  near  and  intermediate  zones.  Due  to  SAM  to  OAM  transformation 
in  the  case  of  a  Rayleigh  scatterer,  the  Poynting  vector  propagates  along  conical  and 
spiral  trajectories  (Fig.  33).  For  larger  spheres,  the  more  complicated  process  of 
scattering  results  in  a  complex,  sometimes  winding  trajectory  of  S. 


Figure  33.  Poynting  vector  field  lines  projected  on  the  plane  perpendicular  to 
the  direction  of  excitation  and  containing  the  center  of  the  scattering  sphere. 
Different  size  spheres,  smaller  than  the  wavelength  X  and  larger  than  X  as 
indicated,  exhibit  a  spiraling  of  the  flow  of  energy  in  the  intermediate 
region.  Most  interesting  is  the  winding  of  Poynting  vector  field  lines  in  due 
to  the  complex  process  of  scattering  from  large  spheres. 


For  instance,  Fig.  33  illustrates  the  projection  of  S  on  the  plane  perpendicular  to  the 
direction  of  excitation  and  containing  the  center  of  the  sphere.  Notably,  the  bending  of 
the  Poynting  vector  field  lines  as  seen  in  Fig.  3  determines  large  angles  between  the  S 
and  p,  resulting  in  experimentally  measurable  SOI  effects.  This  intermediate  zone  spans 
up  to  several  X,  and,  therefore,  an  observation  can  be  performed  without  significantly 
disturbing  the  field  near  the  sphere’s  interface.  The  large  intermediate  zone  for  a  sphere  is 
due  to  the  evanescent  fields,  which  in  this  case  exhibit  inverse  power  law  decay,  as 
opposed  to  the  exponential  decay  for  a  planar  interface. 

Observing  the  circulation  of  the  Poynting  vector  in  scattering  from  a  sphere  requires  a 
detection  scheme  sensitive  to  the  local  energy  flow.  One  such  experiment  could  involve 
optical  forces,  where  the  trajectories  of  probe  scatterers  are  analyzed  in  the  proximity  of 
the  sphere.  A  much  simpler  procedure  however  can  be  based  on  a  common  tool  having 
angular  selectivity,  i.e.  a  single-mode  optical  fiber.  The  coupling  efficiency  of  a  single 
mode  fiber  effectively  depends  on  the  local  field  distribution 


(51) 


where  Es  and  H  are  the  scattered  field  and  fiber  mode,  y  is  the  unit  vector  along  the  fiber 
and  r  is  the  electric  field  transmission  coefficient.  The  integration  is  performed  in  the 
plane  o  perpendicular  to  the  fiber. 


Figure  33.  Due  to  spin-orbit  interaction,  a  circularly  polarized  plane  wave 
suffers  a  transversal  shift  A  in  any  plane  X.  This  shift  depends  on  the  direction 
of  observation  0  and  is  a  manifestation  of  the  spin  Hall  effect  in  scattering. 
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Figure  34.  Schematic  of  the  experimental  setup  using  a  single  mode 
optical  fiber  as  a  means  for  sensing  local  power  flow  direction. 


We  used  an  experimental  setup  based  on  an  array  of  single-mode  optical  fibers  (a 
coherent  fiber  bundle)  as  depicted  in  Fig.  34.  A  series  of  measurements  were  performed 
as  the  fibers  were  scanned  in  the  intermediate  region  of  dielectric  spheres,  perpendicular 
to  the  direction  of  propagation  of  the  incident  field.  To  avoid  interferences  with  the 
substrate,  the  spheres  were  suspended  in  oil  matching  the  index  of  the  microscope  slide. 
To  reveal  the  asymmetry  in  the  field  distribution  resulting  from  excitation  with  different 
spins,  we  performed  sequential  measurements  and  recorded  the  corresponding  scattered 
intensities.  Subsequently,  these  values  were  subtracted  to  reveal  the  shift  in  the  direction 
of  S. 

We  have  also  calculated  analytically  the  intensity  coupled  through  single-mode  optical 
fibers  after  scattering  from  a  sphere  and  compared  them  to  the  experimental  results.  The 
remarkable  agreement  is  illustrated  in  Fig.  35  and  demonstrates  not  only  the  asymmetry 
of  the  scattered  fields,  but  also  our  capability  to  infer  a  shift  in  the  perceived  location  of 
the  interaction  volume  (sphere’s  location).  This  constitutes  the  first  observation  of 
transversal  spin  transport  (a  manifestation  of  the  spin  Hall  effect)  in  scattering  from  finite 
size  objects,  using  an  object  rather  than  the  exciting  field  to  localize  the  interaction. 

In  conclusion,  due  to  spin-orbit  interaction,  an  encounter  with  a  gradient  of  refractive 
index  leads  to  transport  of  spin  similar  to  the  electronic  spin  Hall  effect.  Our  experiments 
have  demonstrated  that  transversal  spin  transport  is  possible  even  when  the  symmetry  of 


optical  interaction  is  of  higher  dimensionality.  For  a  wave  in  a  pure  state  of  polarization, 
the  spin-orbit  interaction  results  in  a  spiraling  power  flow  that  is  determined  by  the  extent 
of  the  interaction.  These  results  constitute  the  first  demonstration  that  spin  transport  can 
be  manipulated  and  enhanced  in  confined  geometries,  and  open  the  possibility  for 
developing  new  functionalities  for  sensing  devices. 
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Translation  of  optical  fiber  (pm) 


Translation  of  optical  fiber  (pm) 


Figure  35  Comparison  between  the  analytical  prediction  (a)  and  the  experimental 
results  (b)  for  the  coupled  power  through  a  single  mode  fiber  scanned  across  a 
polystyrene  sphere  of  4.6pm  diameter.  The  graphs  depict  the  difference  between  the 
scattered  intensities  corresponding  to  pure  states  of  excitation  plotted  as  a  function  of 
fiber’s  position  with  respect  to  the  center  of  the  sphere. 


The  significance  of  our  demonstration  is  manifold.  First,  optical  experiments  are  suitable 
tools  to  model  spin  phenomena  that  in  electronic  conduction  may  be  difficult  or  even 
impossible  to  approach.  This  could  lead  to  the  discovery  of  novel  manifestations  of  spin 
transport  in  confined  geometries  where  effects  such  as  loss  of  coherence  and  dissipation 
are  expected  to  be  minimized.  Second,  understanding  subtle  aspects  of  conservation  laws 
in  optical  scattering  should  provide  insights  into  phenomena  such  as  spin  transfer  and 
power  flow  which,  in  turn,  are  essential  for  developing  new  sensing  approaches  at 
nanoscales. 


Manipulating  the  polarization  properties  of  electromagnetic  fields  may  also  have 
consequences  for  controlling  the  subwavelength  behavior  of  optical  forces.  Lastly, 
circular  polarization  is  rather  exotic  in  nature  and  also  robust  in  propagation.  The  ability 
to  distinguish  between  left  and  right  may  provide  unique  communication  channels. 


7.  Conservative  and  nonconservative  optical  torques 

A  fundamental  consequence  of  an  applied  force  is  the  ability  to  induce  torque  with 
respect  to  some  reference  point.  Torques  can  also  be  induced  by  optical  fields.  Several 
concepts  for  optical  spin  motors  or  “ nano-rotators ”  have  been  discussed  based  on  optical 
traps  created  with  circularly  polarized  light  or  vortex  beams  and  relaying  on  object’s 
asymmetry,  absorption,  or  birefringence.  Another  notable  proposal  uses  the  subtle 
interplay  between  conservative  and  nonconservative  forces  in  an  optical  trap  to  create  a 
“nano -fountain”  with  constant  circulation  of  trapped  particles.  The  idea  of  mechanical 
action  of  light  has  been  pursued  for  hundreds  of  years.  The  ability  to  trap  and  maneuver 
small  objects  such  as  microparticles,  polymer  chains,  cells,  etc.  is  undoubtedly  one  of  the 
most  exciting  use  of  what  is  now  commonly  referred  to  as  optical  tweezers.  A  host  of 
applications  are  being  pursued  where  optical  forces  are  employed  for  manipulation, 
measurements,  or  for  creating  and  controlling  new  states  of  matter. 

Two  particles  excited  by  a  common  field  can  form  a  bound  “optical  dimer”  when  they 
arrange  themselves  to  a  stable  position  where  the  radial  forces  acting  on  them  are  zero. 
We  studied  the  canonical  case  of  two  lossless  spheres  that  are  electromagnetically 
coupled.  Electromagnetic  fields  can  induce  conservative  forces  resulting  from  field 
gradients  as  well  as  nonconservative  forces  appearing  due  to  radiation  pressure  and 
gradients  of  phase.  The  calculation  of  optical  forces  acting  on  matter  is  believed  to  be  a 
well-established  routine.  The  same  is  true  for  the  optical  forces  arising  in  basic  OB 
situations,  even  though  care  must  be  taken  in  describing  the  particle-particle  interaction. 


Figure  36  Optical  binding  in  elliptically  polarized  light  E/.  Apart  from  the 
binding  force  FR,  interacting  particles  experience  tangential  forces  Ft-.  Note  the 
existence  of  differential  forces  AF  leading  to  individual  spinning  in  addition  to 
common  orbiting  of  particles  around  the  system’s  center  of  mass. 


Let  us  examine  the  system  of  two  identical  spherical  particles  illuminated  by  a  plane 
wave  propagating  perpendicularly  to  the  radius -vector  connecting  the  centers  of  the 
particles,  as  shown  in  Figure  36.  Forces  are  generated  on  the  spheres  due  to  the  three- 
dimensional,  polarized  field  established  as  a  result  of  scattering.  Due  to  symmetry  in  the 
x-y  plane,  the  force  acting  on  each  particle  can  be  decomposed  into  radial  (binding)  and 
tangential  (rotational)  components.  There  is  also  a  scattering  force  along  k,  but  its  effect 
is  identical  for  the  two  particles  and  does  not  hamper  their  transversal  movement. 


In  the  simplest  case  of  bound  Rayleigh  particles  (dipoles),  the  force  acting  on  one 
scatterer  can  be  estimated  as 
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where  a  is  the  scatterer’ s  polarizability,  E  is  the  electric  field,  u=x,y,z,  and  *  denotes 
complex  conjugate.  The  field  E  is  found  by  solving  self-consistently  the  system  of 
equations  that  takes  into  account  the  mutual  interaction  between  the  particles 

E(r, )  =  E/(r,)  +  GaE(r2 ), 

E(r2 )  =  E;  (r2 )  +  GaE(rj ). 


(53) 


In  Eq.  (53),  r,  and  r2  represent  the  dipole  locations,  E,  (r)  is  the  incident  field,  and  the 
tensor  G(lr,-r2l)  is  the  inter-dipole  propagator.  The  field  derivative  is  then  calculated 
to  obtain  the  final  expression  for  the  force  in  Eq.(52). 


A  popular  way  for  evaluating  the  derivative  dE/ du  is  to  differentiate  the  final  solution  of 
the  system  of  equations  (53).  By  doing  so,  however,  the  results  contradict  the  calculation 
of  time-averaged  forces  based  on  the  well-established  formalism  of  momentum  flux 
tensor  (Maxwell  stress  tensor).  The  correct  way  of  evaluating  the  derivative  dE/du  is  to 
differentiate  with  respect  to  either  r,  or  r2  directly  in  Eq.  (53).  Nevertheless,  it  is 
interesting  to  note  that  the  way  dE/du  is  calculated  has  a  minor  effect  for  the  radial, 
binding  force;  this  is  perhaps  the  reason  this  inconsistency  has  not  been  noticed  before. 
When  evaluating  the  tangential  forces,  however,  there  are  situations  where  the  way  the 
calculation  of  field  derivatives  is  conducted  becomes  important  as  it  will  be  demonstrated 
here. 


Using  Eqs.  (52)  and  (53)  one  can  now  evaluate  the  radial  and  tangential  forces  to  be 
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where /l  =  2exp(ikR)(-ikR  +  l)/ R3 ,  and  //  =  exp ( ikR)(k 2 R 2  +ikR-l)/ R3  are  eigenvalues 
of  G,  R  is  the  distance  between  particles,  and  k  is  the  wavenumber.  E[ ,  F,'  are  the 

components  of  incident  field  perpendicular  and  parallel  to  the  separation  vector.  We  can 
now  proceed  to  examine  the  effect  of  the  incident  polarization. 


7.1  Optical  binding  with  linearly  polarized  light 

Because  the  depth  of  the  potential  wells  in  the  stationary  points  depends  on  the  incident 
polarization,  the  system  of  optically  bound  particles  tends  to  orient  itself  such  that  it 
occupies  the  most  energetically  favorable  position.  This  may  be  interpreted  as  the  result 
of  an  effective  anisotropy  of  the  system  ofspheres  bound  as  a  consequence  of  optical 
interaction.  When  the  interaction  is  weak  (Aa«  1,  jua«  1),  Eq.  (54b)  simplifies  to 

(F7.)=Ici;I2IE/  I2  cos(2#)Re((/t  -ju)/R) ,  where  6  is  the  angle  between  polarization  and 

separation  vectors.  We  note  that  the  tangential  force  varies  in  space  proportionally  to 
cos (kR)  or  sin (kR)  having  the  same  periodicity  as  the  radial  (binding)  force.  The 
tangential  force  acting  on  a  dipole-like  particle  is  zero  when  the  field  polarization  is  along 
or  orthogonal  to  the  separation  vector. 

For  systems  of  larger  particles  there  are  no  simple  analytical  results;  one  has  to  go 
beyond  the  simple  dipole  approximation  and  use  numerical  procedures.  The  method  of 
choice  is  the  Coupled  Dipole  Approximation  (CDA)  where  a  continuous  volume  is 
modeled  as  an  array  of  point  dipoles  responding  to  both  the  excitation  and  mutual 
interactions.  For  modeling  OB  effects,  we  use  an  extension  of  CDA  that  accounts  for  the 
interaction  between  particles  without  discretizing  of  space  between  them.  The 
calculations  yield  the  local  field  distribution  from  which  the  forces  acting  on  each 
individual  dipoles  are  found.  Subsequently,  one  may  readily  find  the  corresponding 
torques  rr  =  ^  Rf  xFr(r/.)  by  summing  over  all  dipoles  in  the  system.  Here  R2  is  the 

component  of  the  radius  vector  perpendicular  to  the  chosen  axis.  In  the  case  of  OB 
spheres,  one  may  identify  torques  resulting  in  two  special  rotational  motions:  (i)  sphere 
orbiting  together  around  their  common  center  of  mass  and  (ii)  spheres  spin  about  their 
own  axis. 
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Figure  37  Torques  in  an  optically  bound  system  of  silica  spheres  in  water 
excited  with  a  field  polarized  linearly  at  an  angle  0  with  respect  to  the  optical 
binding  vector:  (a)  orbital  torque  about  the  system’s  center  of  mass  and  (b)  spin 
torque  of  a  sphere  about  its  own  axis.  Results  are  shown  for  spheres  of  different 

sizes  as  indicated  and  the  torques  are  normalized  to  f0  =|E/|2<r/6/23 . 


Typical  results  of  CDA  calculations  are  illustrated  in  Fig.  37.  As  can  be  seen,  there  are  no 
torques  when  the  incident  polarization  is  orthogonal  or  along  the  separation  vector. 
However,  torques  arise  at  any  other  angle  resulting  in  orbital  and  spin  motions.  The 
unexpected  appearance  of  spin  torque  is  due  to  a  gradient  in  the  tangential  force  across 
the  diameter  of  the  spheres  as  shown  in  Fig.  36.  Note  that  the  torque  does  not  reach  its 
maximum  for  0=  tc/4  as  may  have  been  expected.  The  reason  is  that  the  separation  vector 
corresponding  to  a  stable  binding  position  is  also  a  function  of  angle  0.  Another 
observation  is  that,  remarkably,  the  spin  torque  can  change  sign  for  different  polarization 
angles.  In  fact,  the  mere  existence  of  these  spin  torques  is  a  significant  result, 
demonstrating  that  OB  interaction  can  lead  to  rotations  of  lossless  dielectric  objects. 

The  torques  in  Fig  37  are  mostly  determined  by  gradient  forces  and,  hence,  determined 
by  the  conservative  part  of  the  total  force.  In  any  system  with  damping,  the  mechanical 
motion  created  by  a  conservative  force  will  eventually  cease.  It  follows  that  the  OB 
particles  will  eventually  align  perpendicularly  to  the  direction  of  polarization  and  that  the 
time  scales  of  their  motion  will  depend  heavily  on  the  specific  constraints  of  the  damping 
mechanism.  In  the  following  we  will  reveal  other  situations  where  the  nonconservative 
forces  are  the  main  cause  for  such  torques,  which  can  be  controlled  by  the  polarization  of 
the  external  electromagnetic  field. 

7.2  Optical  binding  with  circularly  polarized  light 

Recently,  we  demonstrated  that  scattering  of  circularly  polarized  light  from  a  sphere 
generates  a  spiraling  energy  flow  around  it.  This  effect  arises  from  the  conversion  of  spin 
angular  momentum  of  incident  light  into  orbital  angular  momentum  of  scattered  light. 
One  can  envision  that  a  test  object  placed  in  the  vicinity  of  such  a  sphere  will  experience 
the  radiation  pressure  from  the  curved  power  flow,  causing  the  object  to  move  along 
curled  trajectory.  In  reality,  the  situation  is  complicated  by  the  interaction  between  the 
two  bodies  as  was  discussed  before.  Moreover,  together  with  radiation  pressure,  the  field 
gradient  force  and  the  force  due  to  gradient  of  phase  may  play  a  significant  role.  Thus, 
the  real  outcome  can  only  be  found  by  analyzing  self-consistently  the  problem  of 
electromagnetic  interaction. 

Starting  from  Eq.  (3b)  in  the  simple  case  of  small  non-absorbing  dielectric  particles,  the 
tangential  force  can  be  approximated  to  be 

(Ft)  =  ±  I  a  I2  Re(a)  I E,  I2  [6kR(3-k2R2)cos(2kR) 

(55) 

-(9 -l5k2R2  +  k4R4)sin(2kR)]/ 4R\ 


The  sign  is  determined  by  the  polarization’s  handedness.  It  is  worth  noting  that  the  force 
magnitude  changes  as  a  function  of  R  with  twice  the  frequency  compared  to  the  optical 
binding  force  evaluated  from  Eq.  (54a).  Furthermore,  contrary  to  the  case  of  linear 
excitation,  the  potential  landscape  is  now  circularly  symmetric  as  shown  in  the  inset  of 
Fig.  38.  This  means  that  the  tangential  forces  are  completely  nonconservative  and  create 
a  steady-state  orbital  torque  about  the  system’s  center  of  mass.  In  addition  to  this 
continuous  rotation  around  the  common  axis,  the  particles  will  also  exhibit  a  continuous 


rotation  around  their  own  axes  due  to  the  gradient  of  the  nonconservative  tangential  force 
along  the  radial  direction. 

To  estimate  the  torques  acting  on  larger  particles  we  used  again  the  CDA  numerical 
approach.  A  typical  summary  plot  of  the  orbital  torque  for  the  first  two  stable  bound 
positions  is  shown  in  Fig.  38  as  function  of  particle  radius.  Also  shown,  are  the  analytical 
predictions  of  Eq.  (55)  for  Rayleigh  particles,  which  seem  to  make  a  good  description  up 
to  a  radius  of  about  a  ~  /t  / 1 0 .  As  apparent  in  Fig.  38,  an  interesting  effect  occurs  for 
larger  spheres:  the  orbital  rotation  can  change  its  sense  depending  on  the  particle  size. 
This  change  in  the  direction  of  rotation,  not  present  in  the  case  of  small  particles,  can 
happen  even  when  moving  between  the  different  stationary  orbits.  Our  calculations  also 
indicate  that  for  particles  with  a  «  A ,  the  radial  and  tangential  forces  have  now  similar 
periodicities  as  a  function  of  R  and,  moreover,  the  zeros  of  radial  force  (stationary  orbits) 
and  the  zeros  of  the  tangential  force  occur  approximately  in  the  same  place.  Thus,  a 
slight  modification  in  the  radial  position  of  spheres  can  change  the  direction  of  rotation. 
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Figure  18.  Magnitude  of  orbital  torque  as  a  function  of  the  radius  of  interacting  spheres  for 
the  first  (blue)  and  second  (red)  stationary  orbits.  The  plus  symbols  indicate  regions  where 
the  torque  has  opposite  sign.  The  dashed  lines  indicate  the  analytical  predictions  based  on 
Eq.  (55)  for  Rayleigh  particles.  The  calculations  are  for  silica  spheres  in  water  excited  with 

a  plane  wave  of  intensity  10 mW / fim  .  The  black  line  shows  the  magnitude  of  torque  due 
to  Brownian  force  at  290K.  The  inset  depicts  the  symmetric  potential  energy  landscape  and 
the  trajectory  of  a  bound  particle  due  to  nonconservative  orbital  torques. 

In  addition  to  electromagnetic  interaction,  OB  systems  can  also  be  subject  to 
environmental  effects.  In  fluids  for  example,  the  tangential  force  is  compensated  by  the 
viscous  Stokes  force  which,  for  spherical  particles,  is  Fs  =  67rr/av ,  where  77  and  v  are 
the  dynamic  viscosity  and  the  particle  velocity,  respectively.  For  particles  of  radius 
a  =  0.41/Lm ,  for  instance,  the  values  for  orbital  torques  in  Figure  18  can  be  used  to  show 


that  maximal  angular  velocity  is  Q0  =  2 v/R  «  3.7  rad/s  in  the  first  stationary  orbit.  This 
angular  velocity  is  easily  detectable  experimentally. 

Another  experimental  aspect  is  Brownian  motion.  Directional  motion  due  to  optical 
forces  will  be  affected  by  the  additional  chaotic  movement  associated  with  some  random 

force  =  127rr/akBT .  The  torque  resulting  from  the  Brownian  force  provides  a  useful 

reference  for  the  magnitude  of  orbital  torques.  In  Figure  18,  one  can  clearly  see  that  for 
a  >  0.32m  and  an  optical  intensity  of  10mW/jum2  ,  the  optically  induced  torques 
dominate. 

Due  to  the  complex  interaction,  the  OB  particles  are  also  subject  to  spin  torque  with 
respect  to  the  individual  axes  as  shown  in  Figure  39.  As  can  be  seen,  for  the  chosen 
parameters,  the  spin  torque  increases  with  the  particle  size  but,  similar  to  the  orbital 
torque,  the  sense  of  rotation  is  not  always  the  same.  Examining  the  two  types  of  torques 
in  Figs.  3  and  4,  one  can  see  that  the  spin  and  orbital  torques  have  the  same  direction  for 
small  particles  but  their  behavior  becomes  more  complicated  when  the  sizes  increase. 
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Figure  39  Magnitude  of  spin  torque  Ts  as  a  function  of  the  radius  of  interacting 

spheres  for  the  first  (blue)  and  second  (red)  stationary  orbits.  The  plus  symbols 
indicate  regions  where  the  torque  has  opposite  sign.  The  calculations  are  for  silica 

spheres  in  water  excited  with  a  plane  wave  of  intensity  10 mW  / /urn2  .  The  black 
line  shows  the  magnitude  of  absorption-induced  spin  torque  of  one  silica  sphere 
with  refractive  index  nt  =  1 .59  + 1CT7  i . 

Circular  polarization  can  induce  torques  on  a  small  object  due  to  asymmetry,  absorption 
or  birefringence.  Therefore,  it  is  instructive  to  compare  the  magnitude  of  OB  spin  torque 
with  the  optical  torque  exerted  on  a  particle  due  to  its  intrinsic  absorption.  The  later  can 

be  estimated  as  Tabs  =  \Ej\2  a2  QabjAk  [*],  where  Qabs  is  the  absorption  coefficient. 


Estimations  based  on  typical  values  for  absorption  in  silica  are  shown  in  Fig.  4,  and,  as 
can  be  seen,  spin  torque  dominates  for  the  entire  range  of  particle  sizes.  Notably,  because 
the  OB  spin  torque  does  not  necessarily  have  the  same  direction  as  the  excitation 
handedness  while  the  torque  due  to  absorption  is  always  in  the  same  direction,  the  two 
torques  can  combine  to  increase  or  to  cancel  the  net  rotation.  Our  calculations  show  that 
for  these  two  types  of  torques  to  have  similar  magnitudes,  the  imaginary  part  of  the 
refractive  index  should  be  on  the  order  of  10“3.  Finally,  one  can  estimate  the  typical 
angular  speeds  Q  s  of  spinning  using  the  fact  that  the  torque  imparted  to  a  rotating  sphere 

by  surrounding  viscous  liquid  is  defined  as  Tv  =  -8/r//fl3Qv .  For  a  sphere  of  radius 
a  =  0.41/tm  and  an  intensity  of  \OmW/ jum2  one  finds  that  the  corresponding  spin  angular 
velocity  in  water  is  Qs  =  rj&7rr/a 3  «  32  rad /s . 

These  results  constitute  the  first  demonstration  that  optical  interaction  forces  can  lead  not 
only  to  binding,  but  also  to  complex  rotations.  The  interplay  between  conservative  and 
nonconservative  forces  constitutes  a  new  mechanism  to  induce  torques  on  spherically 
symmetric,  optically  isotropic,  and  lossless  objects. 

We  have  found  that  when  the  incident  field  is  linearly  polarized,  the  torques  are  mostly 
conservative  and  affect  only  the  transient  behaviors.  For  circular  polarization  on  the  other 
hand,  the  nonconservative  torques  are  significant  and  lead  to  nontrivial  phenomena.  In 
particular,  we  have  shown  that  bound  systems  can  rotate  not  only  around  the  common 
center  of  but  also  around  their  own  axes.  In  the  intermediate  case  of  elliptically  polarized 
light,  the  conservative  torque  will  determine  a  transient  orbital  motion,  whereas  the 
nonconservative  one  will  lead  to  a  continuous  spin  rotation.  The  whole  system  can  be 
seen  as  a  ‘nano-mixer’  with  complex  mutual  rotations  of  constituents.  The  direction  and 
speed  of  these  rotations  can  be  dynamically  controlled  through  the  intensity,  state  of 
polarization,  and  spatial  profile  of  the  incident  radiation.  Our  estimations  indicate  that 
effects  are  easily  observable  under  reasonable  environmental  conditions.  Finally,  the 
bound  system  discussed  here  constitutes  a  new  kind  of  “optical  matter”  having  its 
mechanical  properties  strongly  coupled  to  the  exciting  radiation. 
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